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NOTATION AND UNITS 


Notation 

A brief list of symbols is given here. Unlisted symbols are defined in the chapter where they appear. 


A x , Ay, A z 

B T 


D, 


F (F F F ) 
c’ v cx’ cy’ cz J 

K 


P W, (P x . Py, P z ) 


m 


s 

t 


oscillation amplitude along coordinate axis 

defined by Equation 2.44 

defined by Equation 2.45 

defined by Equation 2.46 

control acceleration 

acceleration magnitude 

perturbing acceleration 

orbital eccentricity 

(1) mass 

(2) in Section C.l of Chapter III, the ratio of the mean motions of the Earth and the 
Moon (m = 0.07480133) 

(3) in Chapter VI, the fraction of the total time that the control is on (See page 76.) 
Laplace transform variaoie 

time 


7 l distance ratio (See Figure 2.3 and page 10.) 
p normalized mass parameter (See pages 5 and 6.) 
v see Equation 2.34 

(1) for eccentricity contribution, see Equation 2.6 

(2) for Sun’s effect in the Earth-Moon system, see Equation 3.32 

p see Equation 2.33 

(1) for eccentricity contribution, see Equation 2.5 

(2) for Sun’s effect in the Earth-Moon system, see Equation 3.31 

c o frequency 


vi 



0 


Sun 


© Earth 
3> Moon 

Note: A numerical subscript is added to the symbols B L , C L , D h , and y L for a particular libration point; 
e.g., B L1 refers to L^. 


Units 

Unless other units are specified, the normalized units of the restricted three-body problem will prevail. 
These normalized units are defined on pages 5 and 6. Conversion factors for the normalized units are 
given in Tables 2.7, 2.8, 2.9, and 2.10. (See explanation on page 23 and following pages.) 

Accelerations are frequently expressed in Earth gravity units (g = 9.81 m/sec 2 ). Note that 1 fps/yr - 

1 x 10 _9 g. 




THE CONTROL AND USE OF LIBRATION-POINT SATELLITES* 


by 

Robert W. Farquhar 
Goddard Space Flight Center 


CHAPTER I 
INTRODUCTION 


A. Problem Statement 

The five equilibrium solutions of the restricted three-body problem have intrigued mathematicians for 
many years. More recently these solutions have aroused the interest of engineers as well. The locations of 
these “libration points” relative to a rotating two-body system are essentially constant; three of them are 
located along the line joining the two primary bodies, while the other two form equilateral triangles with the 
two bodies. A satellite placed at one of these points with the proper velocity will be in equilibrium because 
the gravitational and centripetal accelerations acting on the satellite will cancel each other out. However, 
the three collinear points are unstable, and the equilateral-triangle points are only quasi-stable. Therefore, 
some form of translation control will usually be required to maintain the satellite’s position in the vicinity 
of a libration point. This “station-keeping” problem is complicated by perturbative accelerations and addi- 
tional control functions that are dictated by mission constraints. The primary objective of this research is 
to develop general analytical relationships for these control requirements. 

A secondary objective is to examine critically, and in some cases to modify, existing proposals for 
libration-point satellite missions. This aim will be supplemented by new suggestions for possible applica- 
tions of libration-point satellites. These mission analyses are needed to define various control problems 
and to provide motivation for the present research. 

Attitude control, analysis of transfer trajectories, propulsion systems, and measurement techniques are 
beyond the scope of this study and are only discussed in a cursory manner. 


*The information presented herein was submitted as a thesis in partial fulfillment of the requirements for the 
degree of Doctor of Philosophy in the Astronautical Sciences, Stanford University, Stanford, California, July 1968. 
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B. Previous Contributions 


Since their discovery by Lagrange in 1772, the libration points of the restricted three-body problem have 
been a favorite topic for researchers in celestial mechanics. Although this subject has been treated exten- 
sively, with significant contributions by Hill, Darwin, Brown, Moulton, Strdmgren, and others, new papers 
are being published at an increasing rate. With the exception of a few perturbation studies, most of this 
work has little direct bearing on the main problems of the present study. However, some useful results can 
be gleaned from the theoretical literature. An exhaustive list of references can be found in the recent trea- 
tise by Szebehely (Reference 1). A less extensive list, but one that is primarily concerned with libration 
points, is included in the survey paper by Steg and De Vries (Reference 2). 

In sharp contrast to the vast literature concerning the classical problem, only six references that treat 
the problem of control of a libration-point satellite could be found. For an unstable collinear libration point 
of the Earth-Moon system, a pioneering paper by Colombo (Reference 3) has shown that 

(1) To the first order, this point is an exact solution of a restricted four-body problem (Sun, Earth, 
Moon, and satellite). 

(2) It is possible to use a simple linear feedback control for position stabilization. 

(3) It is possible to obtain the acceleration needed for position and attitude stabilization of a libration- 
point satellite by simply varying the magnitude and direction of a solar sail. 

It should be noted that the analysis and the design of the feedback control mentioned by Colombo are not 
treated in his paper. 

An interesting paper by Dusek (Reference 4) has shown that artificial libration points can be generated 
by constant low-thrust forces. Dusek also proved that, by varying these forces as a function of the satel- 
lite’s position, the motion in the vicinity of a point is bounded. A special class of these artificial libration 
points, the “isosceles-triangle points,” will be considered in the present study. 

Dual-axis feedback controls have been treated by Fleming (Reference 5) and by Paul and Shapiro (Ref- 
erence 6). Fleming used an original method to design linear and nonlinear controls for a satellite in the 
vicinity of an equilateral-triangle point. Paul and Shapiro determined general Routh stability conditions for 
all five of the classical libration points. However, the analysis of Paul and Shapiro was limited to three 
rather restricted controls.. 

For a satellite in the vicinity of a collinear libration point, this present author (Reference 7) has shown 
that a single-axis linear control using only range and range-rate measurements will guarantee asymptotic 
stability. The station-keeping cost for this control was given as a function of the measurement noise. Some 
special controls that are needed for an application were also presented. 

In a recent paper, Kononenko (Reference 8) used Pontryagin’s maximum principle to determine the opti- 
mal trajectory between a libration point and some initial point in the vicinity of the libration point. It was 
assumed that the satellite was equipped with an engine of limited power. 

Applications for libration-point satellites have been suggested in a number of papers (References 7 and 
9 to 18). A discussion of some of these proposals will be presented in Chapter VIII. 
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C. Contributions of This Research 

This study is a basic introduction to the libration-point satellite control problem. Most results are pre- 
sented in a general analytical form and can readily be adapted to a particular application. Because the most 
interesting libration-point satellite applications make use of the collinear points, the major portion of this 
research deals with these points. Some preliminary results of this study have already been published by the 
present author (References 7 and 16), but more details are given here. 

The following paragraphs present in summary form the principal results of this research. 

In Chapter II, general equations of motion in the vicinity of a libration point (classical and artificial) 
are formulated, and numerical data for many of the collinear libration points of the Solar System are pre- 
sented to facilitate the application of general results to specific cases. 

Chapter III develops the concept of nominal path control. Also, analytical estimates for corrections to 
a small quasi-periodic orbit around a collinear point are given; effects of nonlinearity, eccentricity, and 
perturbations are considered; and examples of gravitational perturbations, including a new “equilibrium 
solution” of the restricted four-body problem (Sun, Earth, Moon, and satellite), are presented. 

Simple linear feedback controls for collinear and equilateral-triangle points are analyzed, and Floquet 
stability investigations for certain cases are presented in Chapter IV. 

In Chapter V, analytical estimates of station-keeping costs at collinear points are presented, costs are 
given as functions of the measurement noise, and a solar sail control technique at the Earth-Moon collinear 
points is discussed. 

A limit-cycle analysis for an on-off control system at a collinear point is presented in Chapter VI. A 
closed-form solution is obtained for an important special case, and the stability of the limit cycle is also 
examined. An approximate design technique for the general case is presented as well. 

In Chapter VII, a method is presented for stabilizing the position of the mass center of a cable-connected 
satellite at an unstable collinear libration point by simply changing the length of the cable with an internal 
device. 

In Chapter VIII, a novel method for maintaining a continuous communication link between the Earth 
and the far side of the Moon is presented, and other new applications for libration-point satellites are also 
proposed. 
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CHAPTER II 


PRELIMINARY DEVELOPMENTS 


In this chapter, the equations of motion in the vicinity of libration points are formulated. The libration- 
point concept is traced to its origin in the restricted three-body problem, and various assumptions, approxi- 
mations, and modifications are discussed. Relevant numerical data for many of the collinear points are 
listed in the final section. 

A. The Restricted Three-Body Problem 

1. Equations of Motion 

Consider the motion of an infinitesimal body in the gravitational field of two finite bodies (e.g., 
Earth and Moon) revolving around their common bary center in elliptical orbits. It is assumed that the infini- 
tesimal body does not influence the motion of the two larger bodies, hence the terminology, restricted three- 
body problem. The geometry for this problem is shown in Figure 2.1. Here, the two finite bodies are located 
at P 1 and P 2 , while the infinitesimal body is at P 3 . An inertial reference system (£ 0 , Cq) is shown with 
its origin at the barycenter of the two finite bodies, and its ^Q-axis normal to their orbital plane. A rotating 
coordinate system (£, 77, 4), with its £-axis passing through P 1 and P 2 , is also introduced. From the figure, 
it is obvious that 


and 


£ 0 = £ cos 0 - r] sin 0, 

r] Q - ^ sin 6+1 7 cos 6, (2.1) 

<w- 


Differentiation of Equation 2.1 with respect to time gives 

£ 0 = (£ - Orj) cos e - (77 + dg) sin 0 , 

77 q - ~ Orj) sin 0 + (77 + 0 £) cos 0 , (2.2) 

and 

4 = & 

To facilitate analytical manipulations, the following quantities are defined to be unity: 

(1) The mean distance between P 1 and P g . 


5 




Figure 2.1-Geometry for the restricted three-body 
problem. 


and 


(2) The mean angular rate of the finite 
bodies about the barycenter. 

(3) The sum of the masses of the two finite 
bodies. The larger mass at P 1 is taken as 1 - jn, 
while the smaller mass at P g is /x. Therefore, be- 
cause the mass ratio is M R - (1 - ft)/ fx, then 

li= 1/(1 + M b ). 

Unless other units are specified, this normalization 
will be used throughout this research. 

Because the motion of the two finite bodies 
can be obtained from the solution of the two-body 
problem, it is possible to write the distance \P 1 P 2 \ 
and the angular rate 6 as power series in the eccen- 
tricity. These series are (Reference 19) 

\^P 2 \ = ^z-^ sR = 1 + P (2 * 3) 

and 

6 - 1 + v, (2.4) 


where 


p-~e cos (t + (ft) 


+— [1 - cos 2(£ + eft)] + 0(e 3 ) (2.5) 


v = 2e cos (£ + (ft) +-|-e 2 cos 2(£ + <ft) + 0(e 3 ). (2.6) 

The phase angle ft is just the true anomaly at £ = 0. It should be noted that these series converge rapidly 
when e is small, but are divergent when e > 0.6627 .... 


Using Equation 2.2, the kinetic energy of the infinitesimal body can be written as 

T = T [^0 + ’io + £o] = ■ f [<f - + O + + C 2 ] (2.7) 


The potential energy is 


where 

and 

Forming the Lagrangian 


r? = (f- £i) 2+ ?? 2 + C 2 

r l = (<f - <f 2 ) 2 + v 2 + <C 2 - 

L - T - U, the equations of motion are obtained by Lagrange’s method. 


( 2 . 8 ) 

(2.9) 
( 2 . 10 ) 
This gives 
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I - = h + e 2 ( - - f x ) -4^ ' f 2 ). 

r i r 2 


‘rj + 26>£ = -0£ + 6> 2 r? - 


(1 ~ ft) ft' 


r 3 + ' r 3 
1 2 


and 


(1 ~ /*) fi 


... + 

L r i r 2 _ 


c 


( 2 . 11 ) 


For circular orbits (e = 0) of the two finite bodies, Equation 2.11 becomes 

1 l 2 


77 + 2 i ~ 7 )~ 


(1 ~ ft) ft 


1 


and 


4> 


(i - ^) p. 


+ nr\C 


-3 r 3 

L r l r 


( 2 . 12 ) 


2. Jacobian Constant and the Surfaces of Zero Relative Velocity 


It is the purpose of this section to elucidate an important conservation principle of the restricted 
three-body problem. Denoting generalized coordinates and velocities by q’s and <j’s, respectively, the 
kinetic energy of a system can be written as 


where 


and 


Taking the potential energy as 


T(q, q, t) = T 2 +T 1 + T„, 

(2.13) 

n 

T n 5 2_j ai > {9, 

(2.14) 

i, /-I 


n 


t i 5 mi. 

(2.15) 

i=l 


T 0 - T 0 (q, t). 

(2.16) 


U = U(q, t) 


(2.17) 


and using the definitions 
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L = T - U (Lagrangian), 


(2.18) 


and 


E ~ T + U (total energy), 


(2.19) 



(Hamiltonian), 


( 2 . 20 ) 


it follows that 


H = T 2~ T 0 + u = E ~ T 1 - 2Tq. (2.21) 

It is important to note that, in general, the Hamiltonian is not equal to the total energy. Although this result 
is classical, it is often forgotten because “conventional” dynamical systems usually have a kinetic energy 
of the form T = Tg. 

For the restricted three-body problem, it can be seen from Equation 2.7 that the kinetic energy in 
the inertial coordinates is T = T 2 (q), and therefore H - E. If, instead, the rotating coordinates are em- 
ployed, the kinetic energy for a unit mass is T = T g + T ^ + T Q , where 

T 2 =| + i) 8 + £®), (2.22) 

T^k&l-TS, (8-23) 

and 

T 0 =±8 2 ({ 2 + r, 2 )- (2-24) 

Therefore, in the rotating coordinate system, H ^ E. 

If the two finite bodies are in circular orbits, 0-1, and an interesting relation can be found. For 
this special case, the potential energy in the rotating coordinates is U - U(q ) and H = 0. Although the 
Hamiltonian is a constant of the motion, total energy is not conserved* because H i- E. This Hamiltonian is 
related to the “Jacobian constant” C (H = - ^C), which is usually found in the form 

V 2 = 2W - C, (2.25) 

where 

V 2 = + t ) 2 + c 2 

and 

(8.26) 

^ r l f 2 

By determining the Jacobian constant from the initial conditions of the infinitesimal body, bounds can be 
placed on the region of its subsequent motion. Because V 2 is always positive, Equation 2.25 shows that 
motion is possible only so long as 2 W > C. Taking 7 = 0 (zero relative velocity) and making use of Equations 


*This can also be seen by using inertial coordinates, where U - U(q,t). It then follows that H £ 0, and because 
H = E, total energy is not conserved. 
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2.25 and 2.26, surfaces of zero relative velocity corresponding to different values of C can be constructed 
in the configuration space.* 

The curves of zero relative velocity formed by the intersection of these surfaces with the orbital 
plane of the two finite bodies are depicted in Figure 2.2 for the Earth-Moon system = 0.0121507). Here 
the Jacobian constants are C Q = 3.2880, C 1 = 3.1883, C g = 3.1724, C g = 3.0121, C 4 = 2.9980, and C 5 = 
2.9880. As noted earlier, the possible regions of motion for an infinitesimal body with given initial condi- 
tions can be determined from its Jacobian constant. For example, if the body is initially located near the 
surface of the Earth and its Jacobian constant is greater than C Q , it can never leave the region bounded by 
the C Q curve around the Earth. On the other hand, if its Jacobian constant is C 4 , motion can take place 
everywhere except for the region enclosed by the C 4 curve. Figure 2.2 is a contour map of the modified 
potential field of the Earth and Moon in the rotating coordinates. The regions surrounding the Earth and 
Moon in this modified potential field can be viewed as valleys and the point at L 4 as a mountain peak. 
Notice that saddle surfaces exist at L g , and L g . 

Returning to the more realistic case in which the two finite bodies are in elliptical orbits, it can 
be seen that even for rotating coordinates, U = U(q, t) and H f 0. This means that the Jacobian constant 
does not exist in the elliptical case, and the usefulness of Figure 2.2 may be severely limited. However, it 
is the opinion of this author that Figure 2.2 is still useful for intuitive reasoning when'the eccentricity is 
small and time intervals are relatively short. For a discussion of this problem, see References 20 to 22. 

3. Libration Points 

It is well known that there are five equilibrium solutions of Equation 2.11. These equilibrium 
points are located in the orbital plane of the two finite bodies, and their general configuration is shown in 
Figure 2.3. Two of the points, L 4 and L & , form an equilateral triangle with the finite bodies, while the 
remaining three, L v L g , and L g , are collinear. A common name for these points and the one that will be 
used here is “libration points.” 

The existence of these libration points is reasonable from a purely physical -standpoint because 
they are points where the gravitational and centripetal accelerations are balanced 1 in the rotating coordinate 

system. In the circular restficted three-body prob- 
lem, these points will always be stationary. How- 
ever, in the elliptical case, the distances from the 
finite bodies will vary periodically while the 
equilateral-triangle and collinear configurations are 
maintained. A thorough discussion of the equilibrium 
solutions for both the circular and elliptical cases 
can be found in Reference 1. 

From Figure 2.3, it can be seen that the 
collinear libration points are located at distances 
y L R from the finite bodies.! The constant y L is 

♦Details of this construction as well as many fine examples for various values of /x can be found in Reference 1. 
tAn additional subscript is added for a particular point; e.g., y. L1 refers to L*. 



Figure 2.2-Curves of zero relative velocity for the 
Earth-Moon system. 
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obtained by finding the real root of a certain quintic equation. The quintic equations for the three collinear 
points are (Reference 1) 

y\\ ~ (3 ~ pOy£i + (3 - 2n)yf t ^ - ~ n = 0, (2.27) 

yf , 2 + ” ^yL2 + “ 2 ^y S L2 ~ ^L2 “ L2 ~ P = (2.28) 

and 

y| g + (2 + /x)y£ 3 + (1 + 2//)y| 3 - (1 - /z)y| 3 - 2(1 - /z)y L3 - (1 - p) = 0. (2.29) 

Series solutions for y L are also available, but are not recommended for numerical computations at L 1 and 
L g where the convergence is poor. However, abbreviated versions of the series at L 1 and Lg are sometimes 
useful for analytic approximations. The series for L ^ and L g are (Reference 1) 
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+ . . . 


(2.30) 


and 




(2.31) 


The series solution for L g converges very rapidly, and because the real root of Equation 2.29 is near +1, it 
is preferred for numerical computations. This series is given by (Reference 1) 


Yl 3 = 1 + f 



M.2 

84 



3 


761 4 3163 5 

2352 6 7056 6 


30 703 
49 392 



(2.32) 


where e = 7 [i/ 12. 

B. Equations of Motion in the Vicinity of Libration Points 

In this section, the equations of motion in the vicinity of libration points are obtained in a general form 
that is readily adaptable to the perturbation analyses of Chapter III. Although the derivation could proceed 
directly from Equation 2.11, a more basic approach is followed here. All of the classical libration points 
with the exception of L g , which is not too interesting from the standpoint of possible applications, are 
treated. In addition, an artificial libration point that forms an isosceles triangle with the two finite bodies 
is examined. 

1. Collinear Points 

In Figure 2.3, fhe distance R and the angular rate 0 are defined as 

R = 1 + p(t) (2.33) 

and 

<9=1 + v(t), (2.34) 

where the quantities p(t) and v(t) are functions of time that arise from the orbital eccentricity of the finite 
bodies and perturbations from bodies external to the primary two-body system. The unit vectors i and j are 
fixed in the rotating reference frame, which means that 

di N . 

dt ■ (1 + v)> 

and (2.35) 

i = + 

A right-handed triad is formed by taking the cross product k = i x j. 

From Figure 2.3, it is evident that, for a satellite (infinitesimal mass) in the vicinity of either L 1 
or Lg, the position vector relative to the larger body (1 - ^t) is* 

r = [(1 ± y L )(l + p) + xji + yj + zk. (2.36) 

^Whenever double signs appear in this section, the upper sign will hold at L 2 , and the lower sign at L x . 
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(The position vector relative to the libration point is x L - xi + yj + zk.) Differentiation of Equation 2.36 
yields the velocity vector (relative to the larger mass) 


=j^(l ± y L )p + x - (1 + x/)yji +|(1 + z^)j(l ± y L )( 1 + p) + X J + yjj + zk. 


(2.37) 


The position vector for the smaller finite body is 


R = (1 + p)i. 


(2.38) 


For a satellite of unit mass, the kinetic energy is T = • r, and the effective potential energy (for 

motion relative to the larger mass) is given by 


where 


and 


1 R • r 

|r-R| |R[ S _ 

(1 - p) (i _ [(1 ± y L )(! + p) + x ] { 
f i ^\ r 2 (1 + P ) 2 ) 

r i =[(! 1 xlX 1 + p) + x ] 2 + y 2 + z2 

4 ~ [y^ 1 + 1 °) 1 x ] 2 + y 2 + ^ 2 - 



(2.39) 


(2.40) 

(2.41) 


Forming the Lagrangian L = T - U, standard manipulations lead to the equations of motion 
x - 2(1 + v)y + (l ± y L )p - vy = (1 + ^) 2 [(l ± y L )(l + p) + x ] 


a - p) r 


(1 1 } / l)( 1 + p) + x 


T 7 f [^ (1 + p) 1 x ] ■ ( 77 ^ 


y + 2(1 + v)[(l ± y L )p + x ] + v[(l ± y L )(l + p) + x ] 


(l + rfy - 


(1 - /z) 


r 3 t-S 

L r i r 2j 


y. 


and 


z = - 


(1 ~p) P 


-3 r 3 

L r i 2 


(2.42) 


Expanding Equation 2.42 in a power series* including third-order terms, a moderate amount of algebra leads to 

*The complete series converges in a region that is common to the interiors of two spheres. One sphere is centered 
at 1 - p and has a radius of ± 7l)(1 + P)> the other is centered at p and has a radius of \[2y L (l + p)- 
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and 


where 


and 


x - 2(1 + v)y - £(1 + is) 2 + 2(1 + p) 3 # l J x + i>y 

t|c l (1 + p)~ 4 [2x s - (y a + z 8 )] 

+ 2D L (1 + p)“ 5 [2x 2 - 3(y 2 + z 2 j ij x 

- (1 ± y L )h> ~ (1 + p)' 2 [( 1 + p) 2 (l + p) 8 - l]}> 

y + 2(1 + v)x = - t>x +[(1 + is) 2 - (1 + p)~ 3 B Jy ± 3C L (1 + p)~ 4 xy 
-*|-D l (1 + p) _5 [4x 2 - (y 2 + z 2 )j y 

- (1 + y L )[2(l + v)p + (1 + p)*0 , 

z = - (1 + p)~ z B L z ± 3C L (1 + p)“ 4 xz 

- ~d l (i + p)” 5 ^ 2 " (y 2 + z 2 )J z, 


B 


L 


(1 “ p) ^ P 

yi_ 


c 


L 



(1 ~ p) 

( 4 ytf 


D 


L 


(1 " P) , P 

(i * yz .) 6 rl ' 


(2.43a) 

(2.43b) 

(2.43c) 

(2.44) 

(2.45) 

(2.46) 


In Reference 23, it is shown that B L > 1 for all values of p. For p « 1, Equations 2.30 and 2.31 show that 
p = 3 y^. In this instance, the coefficients defined in Equations 2.44 to 2.46 are B L = 4, C L = 3/y L , and 
D L ' 3 ! ?.. 


The accuracy of the derivation can be checked by verifying that the L ^ and L 2 points are equilib- 
rium solutions in the elliptical case. For this to be true, the expressions for p and v given in Equations 2.5 
and 2.6 must vanish when they are substituted into the final terms of Equations 2.43a and 2.43b. Restrict- 
ing this verification to the linear terms, it is only necessary to have 


3p - p + 2v - 0 

and (2.47) 

2p + is - 0. 

It can be seen by inspection that the linear eccentricity terms of Equations 2.5 and 2.6 satisfy Equation 2.47. 

A linearized version of Equation 2.43 is sufficient for many applications. If the quantities p and v 
are also neglected, Equation 2.43 reduces to 
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X - 2 y - ( 2B l + l)x = 0, 


(2.48a) 


y + 2x + (B L - 1 )y = 0, (2.48b) 

and 

z + B l z = 0. (2.48c) 

It is readily seen that the motion perpendicular to the xy-plane is simple harmonic with frequency 
The motion in the xy-plane is coupled, and the characteristic equation is given by 

s 4 - (B l - 2)s 2 - (2 B l + 1 )(B l - 1) = 0. (2.49) 

This equation possesses two real roots, equal in magnitude, but opposite in sign. The two remaining roots 
are pure imaginaries. Because a positive real root exists, the collinear points are unstable. 


2. Equilateral-Triangle Points 


The geometry for satellite motion in the vicinity of the equilateral-triangle points is shown in 
Figure 2.3. Equations 2.33, 2.34, 2.35, and 2.38 are still valid for this section, and the position and 
velocity vectors relative to the larger body are* 


r = 


and 


1 

2 

1 . 


(1 + p) + x| i + 
p + x - (1 + v) 


± ^(1 + p) + y 

± + p) + y 


] + zk 


(2.50) 


j V3. . , N 

1 , x 1 

j± ~2~p + y + (1 + v) 

-p (1 + p) + x 

J) 


j + zk. 


(2.51) 


(The position vector relative to the libration point is r L = xi + yj + zk.) The kinetic energy T =~i • r and 
effective potential energy 


ii B (1+r,)+x ] 


(i + py 


(2.52) 


where 


and 


r i = 


-(1 + p) + x 






+ ±-^-(1 + p) + y + 


-g(i + p) + x 


:J + 1^± ~(1 + p) + yj + 


(2.53) 

(2.54) 


are used to form the Lagrangian L = T - U, which leads to the equations of motion 


^Whenever double signs appear in this section, the upper sign will hold at L 4 , and the lower sign at L 5 . 
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(2.55) 


x - 2(1 + v)y - -~~P ± V3(l + v)p ±~(1 + p)v +■—(! + p)(l + v) 2 


+ (1 + i/) 2 x + vy 


(i - 

~p— i(i + p)+* 

i L 


~~!w ~j (1 + p> +x - 


(i + P y 


. -Vs.. 


. . Vs , 


y + 2(1 + v)x = + -g-p - (1 + i)p — g(l + P)^ ± "g"(l + pXl + 


- i>x + (1 + v) 2 y - 


(1 ~ p) 


4 + jLir-£ 


3 + r 3 

X 2 


X 


(1 + .p) + y 


and 


z - - 


"(1 - p) p 

__d 4 - _ 

r 3 T r 3 
_ 1 r 2_ 


Z. 


A power series expansion* of Equation 2.55 including terms up to the third order is given by 


y = (: 


x - 2(1 + v)y - | (1 + v) 2 -~{1 + p) - 


8 ]x + i 


F ± X?(l - 8,0(1 + P r 3 |y 


3 


2 +-^xy -|§(1 - 2 p)y 2 +4(1 - 2p)z 2 


! ] 


+ ( 1 + p ) ~ 4 [|(l-2^ ) x 
- (1 + p)~ 5 [|f X s ±|§VF(1 - 2 p)x 2 y --^y-xy 2 + |f ^0- “ ^)y 3 

2 1 ^8^ yz2 ] ~\'P ± V3 (1 + v)p ± ^4(1 + p)i> 


+— xz* ± 

O 




[- 


(1 + p) 8 (l + v) 




o-i/q" 

y + 2(1 + x/)x = ±^P(1 - 2/x)( 1 + pF 8 - i> x + ( 

+ J 


2 , 5 
‘ 4 1 


(1 + vY +~(1 + p) 3 y 


- (1 + p)~ 4 [± ^rV3x 2 + ^(! - 2/r)xy ± ~a/ 3 y 2 + ^~z 2 ] 

- (1 + P)" 5 [±||V3(1 “ ^)x 8 - ^x 2 y + ^ VF(1 - 2p)xy 2 

^32 yS ~ 2 P) xz2 “^ yz2 ] T lfP “ C 1 + ^)P 


1 




“■^■(l + p)v + “^(1 + p) 2 |1 ~ (1 + p) d (l + v) 


[■ 


(2.56a) 


(2.56b) 


*The complete series converges in a region that is common to the interiors of two spheres of radii \/2*(l + p). 
These spheres are centered at 1 -/^and/i. 
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and 



(2.56c) 


Some confidence can be put in the accuracy of this result because the eccentricity check of the previous 
section again leads to the conditions given in Equation 2.47. 


If p = v ~ 0 and only linear terms are retained, Equation 2.56 is reduced to 


and 


Q Qv/q 

x-2y-^x+^(l-2p)y = 0, 

y + 2x +-^-(l - 2 n)x -~y = 0, 
z + z = 0. 


(2.57a) 

(2.57b) 

(2.57c) 


The coupling of x and y in Equations 2.57a and 2.57b can be eliminated by a coordinate rotation. This sim- 
plification will be useful when feedback controls are considered in a later chapter. The angle for this rota- 
tion is conveniently determined by expanding the modified potential function to second order (this is suffi- 
cient for linearized equations of motion) and then eliminating the xy-term that causes the coupling. The 
modified potential function to second order is written as 


w =To- 




- V 


=-|(3 - fi) +J- [x 2 ± 2^3(1 - 2 (i)xy + 3y 2 ] -jfz 2 . 


(2.58) 


It is readily verified that the xy-term of Equation 2.58 can be removed by utilizing the coordinate transforma- 
tion (see Figure 2.3)* 


and 


x = x' cos if/ - y' sin if/, 
y = x' sin if/ + y' cos if/, 
Z - z' , 


where 

tan 2tA = + V3 (1 - 2/n). 

With this transformation, Equations 2.58 and 2.57 become 


2W = (3 - [/) + ax' 2 + j3y ' 2 - z' 2 


(2.59) 


(2.60) 


(2.61) 


*The barycenter lies on the x'-axis. 
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and 



x' - 2 y' - ax' = 0, 

(2.62a) 

and 

where 

y 1 + 2x' - /3 y‘ - 0, 
z' y z' - 0, 

of 1 

(2.62b) 

(2.62c) 

and 

a =-g[l + Vl - 3p(l - p)J 

of 1 

(2.63) 


(8 =|[l - Vl - 3,t(l - ft) 

(2.64) 


It can be seen at once that Equation 2.62c is simple harmonic with unit frequency. It is also worth noting 
that for p « 1, a = 3 - 9 p/4 and j8 = 9 p/4. 

The stability of the linearized equations of motion can be examined by writing the characteristic 
equation for the gyroscopically coupled portion of Equation 2.62 

s 4 + s 2 + a/3 = 0 (2.65a) 


or in terms of the mass parameter p 


s 4 + s 2 +-|^p(l - p) = 0. (2.65b) 

The roots of this equation are pure imaginaries if ft < 0.03852 .... and the system is neutrally stable. For 
P > 0.03852 . . . , the roots are complex with one pair possessing positive real parts, and the system is 
unstable. 


3. Isosceles-Triangle Points 


If a satellite moving in the mutual gravitational field of two finite bodies is equipped with a con- 
tinuous thrust device, it is possible to artificially generate new equilibrium positions. This generalized 
libration-point concept was introduced by Dusek (Reference 4), who developed general conditions for the 
existence of these points and examined the stability of a variable thrust satellite in their vicinity. In this 
section, the equations of motion of a satellite at the generalized points that form an isosceles triangle with 
the two finite bodies are derived. This case was not explicitly treated by Dusek. 

The geometry for a general isosceles-triangle point is given in Figure 2.4. Proceeding in the same 
manner as in the two previous sections, the position and velocity vectors are 


and 


r = [(1 + p ) cos 6 + xji + £(1 + p) sin 6 + yjj + zk 

r =<jp cos 6 + x - (1 + */)[(l + p) sin 6 + yjji 

+|p sin 6 + y + (1 + i^)[(l + p) cos 6 + x]|j + zk. 


( 2 . 66 ) 


(2.67) 
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The potential energy is 

(1 - p) 


U 


£(1 + p ) cos 0 + xp 


(l + pY 


( 2 . 68 ) 


where 


and 


Figure 2.4-Isosceles-triangle point. 


r 2 + P) cos ^ + X J 2 

+ [(1 + p ) sin 6 + yj 2 + z 2 

r| = [(1 + p)( cos Q - 1) + xj 2 
+ £(1 + p) sin 6 + yj 2 + z 2 . 


(2.69) 


(2.70) 


This time, the equations of motion are augmented with a planar thrust acceleration (F cx , F cy ) and are given 

by 

x - 2(1 + v)y = - p cos 6 + 2(1 + i/)p sin 6 + (1 + p)i> sin 0 + (1 + p)(l + v) 2 cos 0 

p (1 ~ p) 1 

+ (1 + y) 2 x + vy g— (1 + p) cos 6 + xj 


(1 + p)(cos d - 1) + x - 


+ F , 

r!Y * 


(1 + p) 2 ' * cx 

y + 2(1 + v)x - - p sin # - 2(1 + v)p cos 0 - (1 + p)i> cos (9 + (1 + p)(l + F) 2 sin 0 


(2.71) 


and 


with 


z = - 


- vx + (1 + r/) 2 y 
(1“P) , P 


(1 - p) p 


r 3 + ,.3 
f l f 2 


[(1 + p) sin + yj + F cy) 


r 3 + r 3 
1 l 2 


F = F i + F i. 

c cx cy J 


From Equation 2.71, it can be verified that the isosceles-triangle point is an equilibrium solution if* 


F = F 
c c 


. 6. , 0. 

sin^ji + cos — j 


(2.72) 


(2.73) 


where 


^Equations 2.73 and 2.74 are true so long as p and v are functions of the eccentricity only. In this case, 

P ~ (1 + P) 2 [(1 + ^) 2 (1 + /?) S - l] = 0 and 2(1 + v)p + (1 + p)v =0. If p and v also depend on gravitational perturbations 
from other bodies, these relations no longer hold, and Equations 2.73 and 2.74 will be more complicated. 
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(2.74) 


F = — p — jl — P 2(1 - cos 0)1 3//2 I. 

2(1 + p) 2 (l - cos 0) j L J [ 

In Figure 2.4, it can be seen that this thrust vector is radial with respect to p. Therefore, a simple coordi- 
nate rotation to the (x', y‘ , z ') system will uncouple the thrust acceleration terms in Equation 2.71. Inspec- 
tion of Figure 2.4 shows that this coordinate rotation is 

, e , . o 

x - x' cos—- y' sm — 

and (2.75) 

# . 6 > , e 

y - x' sin^+ y 1 cos^-. 


Taking p-v- 0, and using Equations 2.72 to 2.75, Equation 2.71 is expanded in a power series* in the (x', 
y', z') coordinates. Retaining only linear terms, Equation 2.71 becomes 



x' - 2y' - h^p, d)x‘ -|-[(1 - p) sin 0] y‘ - 0, 

(2.76a) 


y‘ + 2x' - 1-|(1 - p) sin 0jx' - h 2 (p, 6)y' - 0, 

(2.76b) 

and 

z' + h 3 (p, ff)z' = 0, 

(2.76c) 

where 

h^p, 0) =-|(l + cos 0) + 3 cos 0 + 2^2(1 - cos 0)j~ 3//2 |, 

(2.77) 


h 2 (p, 0) =-|(l _ cos 0) - ^<1 - 3 cos 0 - 4^2(1 - cos 0)] _3 ' /2 |, 

(2.78) 

and 

h 3 (p, 0) = 1 + J[2(l - cos 0)J” 3//2 - l|. 

(2.79) 


The coefficients depending on p and 0 in Equation 2.76 are shown as functions of 0 in Figures 2.5 and 2.6. 
Figure 2.5 is drawn for the Earth-Moon system (p = 0.0121507), and Figure 2.6 gives the coefficients for the 
Sun-Earth system ( p ~ 3.0404 x 10~ 6 ). The thrust accelerations in Earth gravity units ( g = 9.81 m/sec 2 ) for 
these two systems are given in Figures 2.7 and 2.8. 

Because the out-of-plane motion is once again simple harmonic, the characteristic equation for the 
coupled portion of Equation 2.76 is used to determine linear stability. This characteristic equation can be 
written in the form 


+ 4 - (h 1 + h 2 ) s 2 + h 1 h 2 ~|(1 - p) 2 sin 2 0 = 0. 


(2.80) 


*The region of convergence for this series is again a region that is common to the interiors of two spheres centered 
at 1 - p and p. The sphere at 1 - p has a radius of and the sphere at p has a radius of 2yjl ~ cos 9. 
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Figure 2.5-Coefficients for an isosceles-triangle point in 
the Earth-Moon system. 


Figure 2.6~Coefficients for an isosceles-triangle 
point in the Sun-Earth system. 



Figure 2.7-Thrust acceleration for an isosceles- Figure 2.8-Thrust acceleration for an isosceles-triangle 
triangle point in the Earth-Moon system. point in the Sun-Earth system. 


For a fixed value of ,i, all of the roots of Equation 2.80 are pure imaginaries when 6 > 0 C , and the system is 
neutrally stable. When d < 0 C , the motion is unstable. The approximate values of 6 C for the Earth-Moon and 
Sun-Earth systems are 38.9° and 2.45°, respectively. 

4. Stability of the Linearized Equations of Motion in the Elliptic Case 

In the three previous sections, the linear stability investigation was limited to the special case, 
p - v - o. For a more realistic situation, where the two finite bodies describe elliptic orbits, the linearized 
equations of motion possess periodic coefficients, and the stability test is more involved. A brief summary 
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of some results for this problem is given below. The results of Bennett (References 24 and 25) are quoted 
here for he has investigated all five of the classical libration points. It should be noted that the stability 
results for the equilateral-triangle points were originally found by Colombo et al. (Reference 26) for the 
Earth-Moon system, and later by Danby (Reference 27) for the general case. 

In his first paper, Bennett (Reference 24) used Floquet theory to examine the stability of the five 
classical libration points. For the collinear points, he found that the degree of instability increased with 
larger eccentricity. His results for the equilateral-triangle points are shown in Figure 2.9. 

Bennett’s second paper (Reference 25) gives an analytical method for finding a characteristic root 
A as a power series in the eccentricity; i.e., A = A Q + A^e + Age 8 +...., where A Q is the root for the circular 
problem. This paper proves that X 1 = 0, which shows that the first correction to a characteristic root is pro- 
portional to the square of the eccentricity. For most of the systems considered in the present research, the 
correction will be very small and will usually be neglected. 

C. Numerical Data for Some of the Collinear Points in the Solar System 

In this section, numerical data for the L 1 and L g points of several Sun-planet and planet- satellite con- 
figurations are presented. These data will be needed to apply general results to certain Solar System points. 
It is important to note that although the derived parameters are given with several significant figures for 
consistency, their accuracy is really limited by the errors in the adopted constants. 

1. Mass-Ratio Dependent Quantities 

The derived constants in Tables 2.1 to 2.6 depend on only one physical parameter, the mass ratio 
of the two finite bodies. The values for the mass ratios listed in column 1 of Tables 2.1 and 2.2 were 

obtained from Reference 28* with the exception of 
the following: 

(1) Mercury: This constant is so poorly 
known that three values are used here to bracket its 
probable range. Recent estimates can be found in 
References 29 to 31. 

(2) Venus: Mariner 2 data (Reference 32). 

(3) Earth: NASA Ad Hoc Standards Com- 
mittee (Reference 33). 

(4) Earth + Moon: Derived from adopted 
values for M/M and M /M . 

vL' vfA j) 

(5) Mars: Mariner 4 data (Reference 34). 

(6) Moon: Mariner 2 data (Reference 32). 

Columns 3 and 4 of Tables 2.1 and 2.2 were obtained 
directly from Equations 2.27 and 2.28, respectively. 

^Original sources are given in this reference. 



Figure 2.9~Stability chart for the equilateral-triangle 
points (Reference 24). Stable region is crosshatched. 
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Table 2.1~Mass parameters and distance ratios for Sun-planet systems. 


Planet 

1 

2 

3 

4 

Mass ratio 

(yy 


y Li 

yL2 

Mercury — 1 

5.5 x 10 6 

1.8182 x 10“ 7 

3.9229 x 10“ 3 

3.9331 x 10“ 3 

Mercury - 2 

6.0 x 10 6 

1.6667 x 10“ 7 

3.8109 x 10“ 3 

3.8206 x 10“ 3 

Mercury - 3 

6.5 x 10 6 

1.5385 x 10“ 7 

3.7106 x 10“ 3 

3.7198 x 10“ 3 

Venus 

4.0859 x 10 5 

2.4474 x 10~ 6 

9.3148 x 10“ 3 

9.3730 x 10” 3 

Earth 

3.32951 x 10 5 

3.0034 x 10~ 6 

9.9704 x 10“ 3 

1.0037 x 10 -2 

Earth + Moon 

3.28906 x 10 5 

3.0404 x 10~ 6 

1.0011 x 10“ 1 2 

1.0078 x 10“ 2 

Mars 

3.098 x 10 6 

3.2279 x 10“ 7 

4.7487 x 10“ 3 

4.7640 x 10" 3 

Jupiter 

1.04736 x 10 3 

9.5387 x 10“ 4 

6.6680 x 10“ 2 

6.9784 x 10“ 2 

Saturn 

3.4997 x 10 3 

2.8566 x 10~ 4 

4.4962 x 10“ 2 

4.6351 x 10“ 2 

Uranus * 

2.29 x 10 4 

4.3666 x 10~ 5 

2.4216 x 10“ 2 

2.4613 x 10~ 2 

Neptune 

1.889 x 10 4 

5.2935 x 10“ 5 

2.5806 x 10“ 2 

2.6258 x 10“ 2 


Table 2.2-Mass parameters and distance ratios for planet-satellite systems. 


Satellite 

1 

2 

3 

4 

Mass ratio 

(y/y 

F 

y lx 

y L 2 

Moon 

81.30* 

1.21507 x 10“ 2 

0.150935 

0.167833 

Io 

3.81 x 10~ 5 

3.8099 x 10” 5 

.023148 

.023511 

Europa 

2.48 x 10~ 5 

2.4799 x 10“ 5 

.020083 

.020355 

Ganymede 

8.17 x 10~ 5 

8.1693 x 10“ 8 

.029781 

.030385 

Callisto 

5.09 x 10" 5 

5.0897 x 10~ 5 

.025474 

.025914 

Mimas 

6.68 x 10“ 8 

6.6800 x 10“ 8 

.0028107 

.0028160 

Enceladus 

1.51 x IO"' 7 

1.5100 x 10“ 7 

.0036877 

.0036967 

Rhea 

4 x 10~ 6 

4.0000 x 10“ 6 

.010966 

.011047 

Titan 

2.48 x 10“ 4 

2.4794 x 10“ 4 

.042920 

.044184 

Triton 

1.28 x 10“ 3 

1.2784 x 10“ 3 

.073339 

.077111 


*y\ 


In Tables 2.3 and 2.4, columns 1 and 2 were computed from Equation 2.44, and columns 3 and 4 
were found by using Equation 2.45. The coefficient D L defined in Equation 2.46 was not included in this 
table because B L and C L are sufficient for most applications. However, some representative values for 
D l are 

(1) Earth + Moon: D L1 = 30240 D Lg = 29243. 

(2) Moon: D L1 = 157.355 D L2 = 91.7003. 
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Table 2. 3 -Con slants for equations of motion for Sun-planet 
systems (see Equation 2.43). 


Planet 

1 

2 

3 

4 

B L1 

B L2 

C L1 

C L2 

Mercury — 1 

4.0237 

3.9766 

766.74 

760.75 

Mercury - 2 

4.0230 

3.9772 

789.23 

783.22 

Mercury — 3 

4.0224 

3.9778 

810.48 

804.48 

Venus 

4.0568 

3.9446 

324.07 

318.07 

Earth 

4.0608 

3.9408 

302.89 

296.89 

Earth + Moon 

4.0611 

3.9405 

301.67 

295.67 

Mars 

4.0287 

3.9716 

633.75 

627.75 

Jupiter 

4.4462 

3.6228 

46.933 

40.984 

Saturn 

4.2905 

3.7412 

68.699 

62.721 

Uranus 

4.1512 

3.8580 

125.88 

119.89 

Neptune 

4.1616 

3.8489 

118.24 

112.25 


Table 2.4-Constants for equations of motion for planet- satellite 
systems (see Equation 2.43). 


Satellite 

1 

2 

3 

4 

b m 

B L2 

C L1 

C L2 

Moon 

5.14760 

3.19042 

21.5117 

15.8451 

Io 

4.1443 

3.8642 

131.59 


Europa 

4.1246 

3.8818 

\ 151.38 


Ganymede 

4.1877 

3.8263 

102.72 

96.734 

Callisto 

4.1594 

3.8508 

119.76 

113.77 

Mimas 

4.0169 

3.9832 

1069.4 

1063.4 

Enceladus 

4.0223 

3.9780 

815.53 

809.52 

Rhea 

4.0670 

3.9349 

275.57 

269.57 

Titan 

4.2764 

3.7525 

71.876 

65.896 

Triton 

4.4960 

3.5873 

42.835 

36.898 


The roots of the characteristic equation at L 1 and Lg, Equation 2.49, are listed in Tables 2.5 and 
2.6. It should be recalled that the two pairs of roots for each point are equal in magnitude, but opposite in 
sign. 

2. Conversion Factors for Normalized Units 

Because virtually all general results will be presented in the normalized units that were introduced 
in Section A.l of this chapter, a set of conversion factors will be needed to apply these results to specific 
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Table 2.5-Roots of characteristic equation for Sun- planet systems 
(see Equation 2.49). 


Planet 

1 

2 

3 

4 

Roots for L 1 

Roots for L g 

Real 

Imaginary 

Real 

Imaginary 

Mercury - 1 

2.5178 

2.0774 

2.4989 

2.0659 

Mercury - 2 

2.5175 

2.0772 

2.4991 

2.0660 

Mercury - 3 

2.5173 

2.0771 

2.4994 

2.0662 

Venus 

2.5309 

2.0854 

2.4860 

2.0580 

Earth 

2.5326 

2.0864 

2.4844 

2.0571 

Earth + Moon 

2.5327 

2.0865 

2.4843 

2.0570 

Mars 

2.5198 

2.0786 

2.4969 

2.0647 

Jupiter 

2.6811 

2.1777 

2.3521 

1.9772 

Saturn 

2.6222 

2.1413 

2.4022 

2.0073 

Uranus 

2.5682 

2.1082 

2.4507 

2.0366 

Neptune 

2.5723 

2.1107 

2.4469 

2.0343 


Table 2.6-Roots of characteristic equation for planet-satellite 
systems (see Equation 2.49). 


Satellite 

1 

2 

3 

4 

Roots for L 1 

Roots for Lg 

Real 

Imaginary 

Real 

Imaginary 

Moon 

2.93206 

2.33439 

2.15867 

1.86265 

Io 

2.5655 

2.1065 

2.4532 

2.0381 

Europa 

2.5577 

2.1018 

2.4604 

2.0425 

Ganymede 

2.5824 

2.1169 

2.4376 

2.0287 

Callisto 

2.5714 

2.1102 

2.4477 

2.0348 

Mimas 

2.5151 

2.0757 

2.5015 

2.0675 

Enceladus 

2.5172 

2.0770 

2.4994 

2.0662 

Rhea 

2.5350 

2.0879 

2.4820 

2.0556 

Titan 

2.6167 

2.1380 

2.4069 

2.0102 

Triton 

2.6997 

2.1892 

2.3368 

1.9681 


systems. The conversion factors for length and angular rate are listed in columns 1 and 3 of Tables 2.7 and 
2.8. Manipulation of these two constants leads to the conversion factors for time, velocity, and acceleration 
given in Tables 2.9 and 2.10. To convert a normalized quantity to ordinary units, multiply by the appropri- 
ate factor. For example, if the normalized acceleration of an artificial satellite in the Earth-Moon system is 
given as F = 1.186 x 10~ 4 , use of Table 2.10 gives 

F - (1.186 x 10" 4 ) x (2.77611 x 10“ 4 g) = 3.292 x 10~ 8 g. 
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Table 2.7-Orbital constants for Sun-planet systems. 


Planet 

la 

lb 

2 

3 

Semimajor axis 
a 

Eccentricity 

e 

Mean motion 
n 

(rad/sec) 

(AU)* 

(km) 

Mercury 

0.387099 

5.7910 x 10 7 

0.20563 

8.2668 x 10" 7 

Venus 

.723332 

1.0821 x 10 8 

.00679 

3.2364 x 10“ 7 

Earth + Moon 

1.000000 

1.4960 x 10 8 

.01673 

1.9910 x 10" 7 

Mars 

1.523691 

2.2794 x 10 8 

.09337 

1.0586 x 10“ 7 

Jupiter 

5.203705 

7.7847 x 10 8 

.04863 

1.6780 x 10“ 8 

Saturn 

9.580337 

1.4332 x 10 9 

.05099 

6.7152 x 10" 9 

Uranus 

19.14103 

2.8635 x 10 9 

.04579 

2.3775 x 10" 9 

Neptune 

30.19825 

4.5176 x 10 9 

.00456 

1.1998 x 10" 9 


*AU = 1.49599 x io 8 km, Mariner 2 data (Reference 32). 


Table 2.8-Orbital constants for planet-satellite systems. 


Satellite 

1 

2 

3 

Mean distance 
from planet 
a 

(km) 

Eccentricity 

e 

Mean motion 
n 

(rad/sec) 

Moon 

3.84405 x 10 5 

0.05490 

2.66170 x 10“ 6 

Io 

4.216 x 10 5 

.0000 

4.1106 x 10" 5 

Europa 

6.708 x 10 5 

.0003 

2.0478 x 10" 5 

Ganymede 

1.070 x 10 6 

.0015 

1.0164 x 10“ 5 

Callisto 

1.882 x 10 6 

.0075 

4.3575 x 10" 6 

Mimas 

1.854 x 10 5 

.0201 

7.7165 x 10" 5 

Enceladus 

2.379 x 10 5 

.00444 

5.3073 x 10“ 5 

Rhea 

5.267 x 10 5 

.00098 

1.6098 x 10" 5 

Titan 

1.221 x 10 6 

.02890 

4.5607 x 10" 6 

Triton 

3.534 x 10 5 

.000 

1.2374 x 10" 5 


Multiplication of the distance ratios in Tables 2.1 and 2.2 by the appropriate length factor fur- 
nishes the values in columns 1 and 2 of Tables 2.11 and 2.12, which give the mean distances between the 
center of the smaller finite body and the collinear points. Comparison of columns 2 and 3 for the planet- 
satellite systems shows that the and points are sometimes located very close to the satellite’s sur- 
face. The variation of the libration-point distance because of the orbital eccentricity is given for a few 
systems in columns 1 to 4 of Table 2.13. 
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Table 2.9— Conversion factors for normalized units for Sun-planet systems. 


Planet 

1 

2 

3 



4 

Time 

(days) 

Velocity 

(m/sec) 

Acceleration 

(m/sec 2 ) 

G 8 ) 

Mercury 

14.001 

4.7872 x 10 4 

3.9575 x 10” 2 

4.0342 x 10" 3 

Venus 

35.762 

3.5021 x 10 4 

1.1334 x 10~ 2 

1.1554 x 10" 3 

Earth + Moon 

58.132 

2.9785 x 10 4 

5.9301 x 10“ 3 

6.0450 x 10" 4 

Mars 

109.34 

2.4130 x 10 4 

2.5543 x 10" 3 

2.6038 x l(T 4 

Jupiter 

689.73 

1.3063 x 10 4 

2.1920 x 10" 4 

2.2345 x 10“ 5 

Saturn 

1723.6 

9.6243 x 10 3 

6.4629 x 10" 5 

6.5881 x 10" 6 

Uranus 

4868.1 

6.8081 x 10 3 

1.6187 x 10" 5 

1,6500 x 10" 7 

Neptune 

9646.7 

5.4202 x 10 3 

6.5031 x 10“ 6 

6.6291 x 10“ 8 


Table 2.10-Conversion factors for normalized units for planet-satellite systems. 


Satellite 

1 

2 

3 

4 

• '.'tV 

I 

Velocity 

(m/sec) 

Acceleration 

(m/sec 2 ) 

(g) 

Moon 

4.34838 

1.02317 x 10 3 

2.72336 x 10" 3 

2.77611 x 10" 4 

Io 

.28157 

1.7330 x 10 4 

7.1237 x 10“ 1 

7.2617 x 10" 2 

Europa 

.56519 

1.3737 x 10 4 

2.8131 x 10" 1 

2.8675 x 10" 2 

Ganymede 

1.1387 

1.0876 x 10 4 

1.1055 x 10" 1 

1.1269 x 10" 2 

Callisto 

2.6561 

8.2008 x 10 3 

3.5735 x 10" 2 

3.6427 x 10" 3 

Mimas 

.14999 

1.4306 x 10 4 

1.1040 

1.1253 x 10" 1 

Enceladus 

.21808 

1.2626 x 10 4 

6.7011 x 10" 1 

6.8309 x 10" 2 

Rhea 

.71898 

8.4787 x 10 3 

1.3649 x 10" 1 

1.3913 x 10" 2 

Titan 

2.5378 

5.5686 x 10 3 

2.5396 x 10" 2 

2.5888 x 10" 3 

Triton 

.93533 

4.3731 x 10 3 

5.4114 x 10" 2 

5.5162 x 10" 3 


3. Comparison of L g Distance and Extent of Planetary Shadow 

It is interesting to compare the distances of the apex of a planet's shadow cone and its U g point. 

The mean distance between a planet and the apex of its geometrical shadow cone is given by 


ds = ( r ®/ r p 


i)' 


(2.81) 


where a is the semimajor axis of the planet’s orbit, R p is the planet’s radius, and Rq- 6.957 x 10 5 km 
(Reference 28). Equation 2.81 was used to compute the planetary shadow distances given in column 4 of 
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Table 2.11-Distance between the smaller finite body and a collinear point 

for Sun-planet systems. 



la 

lb 

2a 



2b 


Mean distance between 

Mean distance between 

L Idillvl 

planet and Lj 

planet and Lg 


(AU) 

(km) 

(AD) 

(km) 

Mercury - 1 

0.001519 

2.2717 x 10 5 

0.001523 

2.2777 x 10 5 

Mercury - 2 

.001475 

2.2069 x 10 5 

.001479 

2.2125 x 10 5 

Mercury - 3 

.001436 

2.1488 x 10 5 

.001440 

2.1541 x 10 5 

Venus 

.006738 

1.0079 x 10 6 

.006780 

1.0142 x 10 6 

Earth 

.009970 

1.4916 x 10 6 

.010037 

1.5015 x 10 6 

Earth + Moon 

.010011 

1.4976 x 10 6 

.010078 

1.5077 x 10 6 

Mars 

.007236 

1.0824 x 10 6 

.007259 

1.0859 x 10 6 

Jupiter 

.34699 

5.1909 x 10 7 

.36314 

5.4325 x 10 7 

Saturn 

.43075 

6.4439 x 10 7 

.44406 

6.6431 x 10 7 

Uranus 

.46352 

6.9342 x 10 7 

.47112 

7.0480 x 10 7 

Neptune 

.77931 

1.1658 x 10 8 

.79295 

1.1863 x 10 8 


Table 2.12-Distance between the smaller finite body and a collinear 
point for planet-satellite systems. 


Satellite 

1 

2 

3 

Mean distance 
between satellite 
and Lj 
(km) 

Mean distance 
between satellite 
and Lg 
(km) 

Radius of 
satellite 
(km) 

Moon 

58020 

64 516 

1738 

Io 

9759 

9912 

1620 

Europa 

13 471 

13654 

1415 

Ganymede 

31866 

32 511 

2450 

Callisto 

47 942 

48 770 

2285 

Mimas 

521 

522 

350? 

Enceladus 

877 

879 

<470 

Rhea 

5776 

5818 

675 

Titan 

52405 

53 949 

2475 

Triton 

25 918 

27 251 

2250 


Table 2.14. The assumed values for the planetary radii listed in column 1 were obtained from Reference 28 
with the exception of the following: 

(1) Mercury: Reference 35. 

(2) Venus: Reference 36. 

(3) Mars: Mariner 4 data (Reference 37). 
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By comparing columns 3 and 4 of Table 2.14*, it can be seen that the Lg point is located beyond the plane- 
tary shadow cone for Mercury, Venus, and Earth. To an observer placed at any of the Lg points for these 
planets, the Sun would always be partially eclipsed. At the L 2 points of Mars, Jupiter, Saturn, Uranus, and 
Neptune, the Sun is totally eclipsed (neglecting atmospheric refraction). A “critical" planetary radius, 
which would make the length of the planetary shadow equal to the L 2 distance, is given in column 2. 


Table 2.13-Distance variation between the smaller finite body and a collinear point. 


Body 

1* 

2* 

3* 


Min. distance 
between body 
and 
(km) 

Max. distance 
between body 
and L ^ 
(km) 

Min. distance 
between body 
and Lg 
(km) 

Max. distance 
between body 
and Lg 
(km) 

Moon 

Mercury - 2 
Earth + Moon 
Mars 
Jupiter 

5.4835 x 10 4 
1.7531 x 10 5 
1.4726 x 10 6 
9.8137 x 10 5 
4.9384 x 10 7 

6.1205 x 10 4 
2.6607 x 10 5 
1.5227 x 10 6 
1.1835 x 10 6 
5.4433 x 10 7 

6.0974 x 10 4 
1.7576 x 10 5 
1.4825 x 10 6 
9.8449 x 10 5 
5.1683 x 10 7 

6.8058 x 10 4 
2.6674 x 10 5 
1.5329 x 10 6 
1.1873 x 10 6 
5.6967 x 10 7 


*d m i n — d(l e )> d max - d(l + e). 


Table 2.14-Lg distance and extent of planetary shadow. 


Planet 

1 

2 

3 

4 

Mean radius 
of planet 
(km) 

Critical radius 
of planet 
(km) 

Mean distance 
between planet 
and Lg 
(km) 

Mean distance 
between planet 
and apex 
of shadow 
(km) 

Mercury - 1 

2480 

2725.6 

2.2777 x 10 5 

2.0717 x 10 5 

Mercury — 2 

2480 

2647.8 

2.2125 x 10 5 

2.0717 x 10 5 

Mercury - 3 

2480 

2578.3 

2.1541 x 10 5 

2.0717 x 10 5 

Venus 

6120* 

6460.2 

1.0142 x 10 6 

9.6036 x 10 5 

Earth 

6371 

6913.4 

1.5015 x 10 6 

1.3826 x 10 6 

Earth + Moon 

6371 

6941.4 

1.5077 x 10 6 

1.3826 x 10 6 

Mars 

3390 

3298.5 

1.0859 x 10 6 

1.1162 x 10 6 

Jupiter 

69750 

45382.0 

5.4325 x 10 7 

8.6745 x 10 7 

Saturn 

58170 

30818.1 

6.6431 x 10 7 

1.3077 x 10 8 

Uranus 

23750 

16712.2 

7.0480 x 10 7 

1.0121 x 10 8 

Neptune 

22400 

17 800.5 

1.1863 x 10 8 

1.5030 x 10 8 


*Radius taken to top of cloud layer. 


*A comparison of the mean distances is sufficient because both distances fluctuate in the same manner for 
elliptical orbits. 
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CHAPTER 111 


PERTURBATIONS AND NOMINAL PATH CONTROL 


The restricted three-body model for the motion of a libration-point satellite is an idealization of the true 
physical situation. A more realistic model must account for the gravitational perturbations of other bodies 
and solar radiation pressure. In this chapter, the effect of these disturbances on a libration-point satellite 
is considered in some detail. Quantitative results are given for a few specific cases, and the concept of 
nominal path control is introduced. 

A. Additional Accelerations Acting on a Libration-Point Satellite 

1. Gravitational Perturbations 

The analysis of the motion of a libration-point satellite is more intricate when gravitational pertur- 
bations of another finite body are introduced. The addition of this fourth body affects the motion of the 
other finite bodies, as well as the satellite, and great care must be taken to use a proper mathematical 
model. A recent discussion of certain models for a “restricted four-body problem” has been given by Mohn 
and Kevorkian (Reference 38). 

When Brown and Shook (Reference 39) calculated the perturbations of a Trojan asteroid* by Saturn, 
they found that the indirect effect produced by Saturn in altering Jupiter’s motion was greater than its direct 
effect. This case is typical of perturbations of this type, and the examples given below show that direct 
and indirect effects are usually comparable. 

2. Solar Radiation Pressure 

The magnitude and direction of the acceleration caused by solar radiation pressure is somewhat 
random for a nonoriented satellite with an odd geometrical shape and variable reflectivity, absorptivity, and 
transmissibility properties. However, for an attitude-controlled satellite with a simple configuration and 
uniform surface properties, the variation of this acceleration is largely predictable. Although this chapter 
views radiation pressure as an undesirable perturbation, it will be seen in Chapters V and VI that this effect 
could be used advantageously for satellite position control. 

The magnitude of the acceleration acting on a plane surface of mass m is given byt 


*The Trojan asteroids are located in the vicinity of the equilateral-triangle points of the Sun-Jupiter system, 
tit is assumed that the incident solar radiation is perpendicular to the plane surface. 
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(3.1) 



where is the solar pressure at 1 AU from the Sun (using the solar constant given in Reference 12, = 

r. p W 

4.50 x 10 newtons/nr), d is the distance from the Sun in AU, A is the area of the plane surface, and C 

P 

is a coefficient that depends on the surface properties. For a perfect absorber (black body), C = 1, while a 
value of C p = 2 corresponds to a perfect reflector. The area-to-mass ratios range from about 0.001 m 2 /kg 
for a very dense satellite to 0.1 m 2 /kg for a satellite with very large solar panels. Assuming an average 
value of (A s /m) = 0.01 m 2 /kg and taking C p = 2, d s = 1 AU, the acceleration from radiation pressure is 
roughly K p = 9.17 x 10~ 9 g. 

3. Thrust Control 

Satellite position control is usually achieved with low-thrust devices. These devices could be 
cold-gas jets, ion thrustors, or some other suitable system. The thrust may be applied either continuously 
or discontinuously, depending on the control requirements. Detailed analyses of the controlled motion of a 
libration-point satellite are presented in the remaining chapters. 

B. Control About Nominal Path 

Although the concept introduced in this section is applicable to all the libration points, the discussion 
here will be limited to the eollinear points. The modification of this technique for the other points is 
straightforward. 

1. Basic Strategy 

The addition of perturbation and control terms to the linearized equations of motion at the eollinear 
points, Equation 2.48, gives 


X - 2y - (2 B l + l)x = P x (t) + F cx , 

(3.2a) 

y + 2x + (B l - 1 )y = P y (t) + F cy , 

(3.2b) 

z + B L Z = P zM + F cz’ 

(3.2c) 


where P(t) is a periodic perturbing acceleration and F c is a control acceleration. Because of the natural 
instability of a satellite in the vicinity of a eollinear libration point, some thrust control for station keeping 
will always be necessary. Additional thrust control could be used to cancel the perturbing acceleration. 
However, it is more economical to control about a suitable periodic perturbed path. This “nominal path” is 
determined by simply solving for the forced response of Equation 3.2 with F c = 0. The forced response to a 
sinusoidal input is well defined, except for exact resonance, even though the uncontrolled system, is 
unstable. Denoting the nominal path by (x n , y a , z n ) and using the coordinate transformation 

£=X-X„, 


30 



(3.3) 


r) = y-y n > 

and 

f= z ~ v 

the equations of motion relative to the nominal path are 

f-2 I j-(2B L + l)f=F cx , 

, + 2f+(B L -l)7) = F cy , (3.4) 

and 

(+B L t=F„. 

which is identical to Equation 3.2 with P(t) absent. This result could have been stated earlier by simply 
using the principle of linear superposition. 

The foregoing development is an approximation that is strictly valid only so long as the nominal 
path is very small. Equations of motion relative to a moderately large nominal path are derived below. 


2. Periodic Orbits 


It is well known that families of unstable periodic orbits can exist in the vicinity of the collinear 
libration points (Reference 1). To satisfy certain mission constraints, it is sometimes necessary to insert 
a satellite into one of these orbits. For moderate amplitude orbits located in the xy -plane, a good approx- 
imation can be obtained from Equations 3.2a and 3.2b with P(t) = F c = 0 . If certain initial conditions are 
satisfied, only the oscillatory mode is excited, and the satellite will follow a periodic orbit about the 
libration point, as shown in Figure 3.1. The equations for this orbit are (Reference 1) 


x„ = 4 X sin a B t 
and 

y n = cos co n t, 


where 


2co, 


A x = 


+ (2B l + 1) 


A y S kA y 


(3.5) 


(3.6) 



Figure 3.1-Periodic orbit around a collinear point. 


and co n is the magnitude of an imaginary root of 
Equation 2.49. It should be noted here that Equation 
3.5 is only a first approximation to the true orbit, 
and higher order corrections are introduced when 
nonlinearities, eccentricity, and perturbations are 
taken into account. The corrected orbit is still 
bounded, but it is only quasi-periodic. 

3. Higher Order Corrections for a Periodic Orbit 

The determination of the corrections to the 
orbit of Equation 3.5 is a tedious and lengthy 
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process. Fortunately, the computation of a few higher order terms is usually sufficient because the remain- 
ing small accelerations do not significantly add to the station-keeping cost. In this section, only second- 
order terms are considered. Although the corrections may be quite small, the station-keeping cost could be 
significantly increased if they are disregarded. 

The method of successive approximations (Reference 40) can be used to determine the second-order 
correction to a periodic orbit around a collinear point. The equations of the corrected orbit to second order 
can be written in the form 


and 


x n = X „1 + x n2 


= y n i + y 


n 2‘ 


(3.7) 


where x fll and y nl are given by Equation 3.5. That is, 


and 


x nl = A xi sin V 


y„ 1 = A yl cos V 


(3.8) 


with A xl = kA yl . Substitution of Equation 3.7 into Equations 2.43a and 2.43b with P (t) included and dele- 
tion of all terms higher than second order yield two coupled linear differential equations for x flg and y ng . 
A yl , p, v, and P(t) are treated as first-order quantities, and coupling between the eccentricity and gravita- 
tional perturbations is neglected. With these restrictions, the second-order effects of nonlinearity, eccen- 
tricity, and perturbations can be treated separately. An analysis of this type is usually adequate for order 
of magnitude estimates of the different effects. 


Frequency corrections are also present and can be computed by expressing the frequency as 

"a = ^nO + "nl + M n2 + ' * * ' 


(3-9) 


Fortunately, a> nl - 0 for the cases considered below, and the frequency is well approximated by co n = (t> nQ . 
Response to a periodic input 

It will be necessary to calculate the forced response to an input of the form 
P(t) = K x cos (cot + 0)i + K y sin ( cot + 0) j . 

The linearized equations of motion near a collinear point with this input are 

x - 2y - (2B l + l)x = K x cos (cot + <fi) 

and (3.10) 

y + 2x + (B l - l)y = K y sin (cot + 0). 

A particular solution of Equation 3.10 can be written as 


and 
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Xp = A x cos (cot + 
y p = A y sin (cot + 0), 


(3.11) 



where 


V 


[-« 8 + (B l - 1 )] 





2 w \-(o 2 - (2 B l + 1 )J 


b 



[® 4 + (B h - Z)o> 2 - (2 B l + 1 )(B l - 

»] 


This result will be used repeatedly below. 


(3.12) 


Nonlinear correction 


Employing the procedure outlined above, it can be deduced that the differential equations for 
the nonlinear correction are 



X „2 - 2 y n 2 - ( 2B L + D x n2 = C x + K x cos 2 <V 


and 

y»2 + 23i „2 + ( b l ~ = K y S™'Zo> n t, 

(3.13) 

where 

c x = ± 2 C ^- k ^) A yv 

(3.14) 



(3.15) 

and 

K y = 4 C L kA2 yv 

(3.16) 

Therefore 

x n2 = x c2 + A x2 cos 


and 

7*2 = A y2 sin 2< V' 

(3.17) 

where A x2 

and A 2 are computed by making use of Equations 3.10 to 3.12, and x cg = -C X /(2B L 

+ 1). As an 

example, consider the L g point of the Earth-Moon system. For this point, B L2 = 3.19042, C L2 

= 15.8451, 


co n - 1.86265, k = 0.343336, and 


C x = 9.0821 A* v 


K y = 8.1603 A 2 V 


x c2 - -1.2305 A* v 


K x = 14.686 A* v A x2 = -0.57443 A 2 V 


A y2 - -0.33203 A yl . 


Choosing ^yl ® .02 = 7688.1 km, these relationships give 


C x = 1.0085 x 10 6 g, x c2 = -189.20 km. 
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K x = 1.6307 x 10" 6 g, A x2 = -88.33 km, 

Ky = 9.0617 x 10 ~ 7 g, A y£ = -51.05 km, 

and A xl = 0.343336 A yl - 2639.6 km. 

It should be mentioned that Plummer (References 41 and 42) has calculated the nonlinear 
corrections to third order for all of the collinear points. Unfortunately, Plummer’s solutions are limited to 
the case where p = 1/11. 

Eccentricity correction 

For the eccentricity correction, the differential equations are 

X n2 ~ ^02 “ (2 B l + X)x„ s = 2vy nl + vy al + (2v - 6 B L p)x nl 
and (3.18) 

y»2 + Si n2 + ( B L - Uy„ 2 = - 2^ nl - ix al + (2* + 3B L p)y nV 

where the expressions for p and v are obtained from Equations 2.5 and 2.6. Taking cf> - 0 and neglecting 
terms of 0(e 2 ), these quantities are given by 

p - - e cos t 

and (3.19) 

v - 2e cos t. 


Using Equations 3.6, 3.8, and 3.19, Equation 3.18 becomes 



*„2 ~ 2y nS “ (2 B l + l)x„ 2 = K s sin (<u n + l)t + K' x sin (<a a - l)t 


and 




(3.20) 


y n2 + 2x n2 + ( B L ~ ^ n2 = K y cos + ^ + K 'y cos ~ 1 ) £ > 


where 

«x = •[- 

2co n - 1 + k(W L + 2)J; 

4 

Vi’ 

(3.21) 


< = - e 

2co n - 1 - k(3B L + 2)ji 

4 

y v 

(3.22) 


K y = ~e 

fi(2o>„ + 1) + (-| B l - gj 

1 A yl , 

(3.23) 

and 


r 




<D 

1 

II 

W*' 

*(2<“ n -l)+(l B L - 2 ] 

i| v 

(3.24) 


The solutions can be written in the form 

x n2 ~ A x2 sin K + + A x2 sin ("a “ 1 ) £ 
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and 


(3.25) 


y n2 = A y2 cos (oi n + 1 )£ + cos (<u n - 1 )t, 

and A x2 , A' x2 , A y2 , and A' y2 are again determined by using Equations 3.10 to 3.12. Returning to the exam- 
ple at the L 2 point of the Earth-Moon system (e = 0.05490), elementary calculations give 

K x - -0.041311 A yl , A x2 = 0.026893 A yV 

K x = 0.068489 A yV A' xZ = -0.030659 A yl , 

K y = -0.24200 A yV A y2 = 0.065947 A yV 

K' y = -0.20430 A yV A y2 = -0.10469 A yV 

and for A y 1 = 0.02 = 7688.1 km, 

K x = -2.2937 x 10~ 7 g, A x2 = 206.76 km, 

K x = 3.8027 x 10 ~ 7 g, A x2 = -235.71 km, 

K y = -1.3436 x I0“ 6 g, A y2 = 507.01 km, 

K y = -1.1343 x 10“ 6 ^ A y 2 = -804.84 km. 

A comparison shows that the amplitudes of the eccentricity corrections are much larger than the amplitudes 
of the nonlinear corrections. It is shown below that these eccentricity corrections also dominate corrections 
caused by the Sun’s gravitational field and radiation pressure. 

Effect of gravitational perturbations and solar radiation pressure 

Because the analysis of gravitational perturbations is slightly different for every libration 
point, a generalized result for periodic orbit corrections cannot be given. However, specific cases are 
readily analyzed, and the effect of the Sun’s perturbation on a periodic orbit around a collinear point of the 
Earth-Moon system is described in section C of this chapter. 

A general result for solar radiation pressure is even less probable because this effect is 
highly dependent on the spacecraft configuration and orientation. Nevertheless, a highly restricted example 
is examined here for illustrative purposes. Consider a satellite near the h 2 point of the Earth-Moon system. 
The major contribution to the satellite’s cross-sectional area is given by a high-gain antenna with a diam- 
eter of 3 m (A s = 7.0686 m 2 ). If the antenna is always directed toward the Moon as shown in Figure 3.2, 
then the angle if/ of the incident solar radiation will vary sinusoidally. Assuming that all of the incident 
radiation is reflected specularly, the instantaneous acceleration acting on the spacecraft is (Reference 43) 

P x = ± cos 2 xfj , (3.26) 
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Figure 3.2-Geometry for radiation pressure example. 


where K p is given by Equation 3.1, with C p = 2. For 
a satellite mass m = 200 kg, the average accelera- 
tion is | PJ = K p / 2 = 1.62 x 10~ 8 g. Accelerations 
of this magnitude may simply be canceled by thrust 
control. 

4. Linearized Equations of Motion Relative 
to a Nominal Path 

A nominal path can usually be represented 
as a series of successive approximations of the form 


and 


k n 1 


y n 


n 1 


+ x nZ + x nS 


+ y n z + y nS 


+ z n2 + z nS 


+ . . . , 

+ . . . , 

+ . . . , 


(3.27) 


where the second subscript denotes the order of the term. The equations of motion relative to the nominal 
path can be found by substituting Equation 3.3 into Equation 2.43. Only terms involving the coordinates 
(<f, V> 0 will be present because Equation 3.27 is a solution of Equation 2.43 with P(t) included. Assuming 
that p and v are negligible when compared to 3C L x nl , etc., the linearized equations of motion are 
approximately* 

£ ~ 2^7 - (2Pj r, + 1)£ = + |^2x nl ^ - y n ^q - + ^ cx . 

ij + - l)p = ± 3 C L jx^^ + + ^ cy> (3.28) 

and 

4 + B L C=±3C L [ X „i<r+z nl £] + F cz . 

When the nominal path is very small, Equation 3.28 reduces to Equation 3.4. 

For a quasi-periodic orbit, the first-order terms for the nominal path are given by Equation 3.5. If 
the satellite is also performing an out-of-plane oscillation, it is obvious from Equation 2.43c that 

z nl = A zl cos (co z t + a Q ) , (3.29) 

where m z - Vb l and a Q is a phase angle. Substitution of these first-order terms into Equation 3.28 gives 
£-2 rj- (2 B l + 1)£= + 3C l | 2^4 x1 sin <o n tj£- ^ yl cos cj 

~~ C0S *]<} + ^CX’ 


*As in Chapter II, the upper sign holds at L s , and the lower sign at L t . 
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V + 2i+ (B l - 1)t] = ± 3C l |[^ x1 sin 
and 

C + 11 sin w n £ j£ + 

C. Examples of Gravitational Perturbations 

1. Solar Perturbation Near the Earth-Moon Collinear Points 

The effect of the Sun’s perturbation on a satellite at a collinear point of the Earth-Moon system 
was first investigated by Colombo (References 3 and 44). He demonstrated that, to first order, the accelera- 
tion produced by the direct effect of the Sun is completely canceled by the indirect effect of the Sun on the 
Moon. A more complete derivation of Colombo’s solution has been given by Nicholson (Reference 45), who 
also examined the forced motion caused by the higher order terms.* It should be emphasized that Colombo's 
solution is valid only at the libration point, and cannot be applied to a satellite following a periodic orbit 
around the point. 

To calculate the effect of the Sun’s perturbation on a periodic orbit, a somewhat different approach 
is advisable. The method presented in this section uses the lunar theory of De Pontecoulant to account for 
the Sun's indirect effect. This method is first shown to be in agreement with Colombo’s result at the libra- 
tion point, and then used to determine the correction to a periodic orbit around a collinear point. Because 
the mean inclination of the Moon’s orbital plane to the ecliptic is only about 5°, and the eccentricities of 
the orbits of the Earth and the Moon are e' - 0.01673 and e = 0.05490, respectively, the analysis given 
below will be limited to zero-eccentricity, coplanar orbits. 

De Pontecoulant’ s expressions for the motion of the Moon 

A complete listing of De Pontecoulant ’s results for the motion of the Moon is given in Refer- 
ence 46, and an abbreviated version can be found in Reference 47. Literally hundreds of terms are involved 
in these expressions, but if it is assumed that e = e' = 0 and three-dimensional effects are neglected, the 
formulas reduce to 

1 2 , 331 4 T 2 , 19 3 , 125 4 1 0 3 4 „ * 

P = ~~6 m + 288 m " C l8 m ] C ° S C0S 4 "s £ 

+ [jjrm +-||m s J-“COs ca s t -'||ra 8 ~-cos 3 co s t + 0(m 5 ) (3.31) 

and 

ril 2 85 3 , 539 4 1 0 , ,201 4 . , 

v- I + l2 m + 36^ ffl J C0S 2(0 s t + cos 4(U s £ 

- | +^-m 8 j 7^7 cos co s t + m 2 ~ cos 3cu s t + 0(m 5 ), (3.32) 


<V]*7 + [_ A yl cos co « £ ]^| + F cy 
[_ A zl COS (co z t + a 0 )]^|+ F cz . 


*The magnitude of the acceleration from the higher order terms is approximately 2 x 10~ 9 g. A nominal path calcu- 
lation is unnecessary here because the magnitude of this acceleration is comparable to random solar radiation pressure 
effects. 
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where m is the ratio of the mean motions of the Earth and the Moon (m = n'/n = 0.07480133), a/a' is the ratio 
of the semimajor axes of the Earth and the Moon [a/a' = 0.002559 ~ 0(m 2 )], co s = 1 - m = 0.92519867, and 
t = 0 corresponds to new Moon. Differentiation of Equations 3.31 and 3.32 gives 


and 


p - ^2ffl s + ^ sin 2 <x> s t 4 —m** sin 4« s t 

- |jfm + -f^ 2 ]-|7 sin co s t +|| sin 3 co s t + 0(m 5 ), 

p = j4m 2 + cos 2co s t + 6m 4 cos 4co s t 

r 15 , 51 2 1 a , , 225 2 a 0 (1 

~ \j£ m + i6 m Ji7 cos "a* + -Ql m # cos *».* + °( m 

fll 2 L 26 3 , 142 4 1 . 0 , 201 4 . . 

v - ~ I — nr + — m + — ^-nr sin 2 co s t sin 4co s t 

, [~15 , 63 2 I a • * 135 2 a o * , r\/ 5\ 

+ [T m+ T ra Ja 7 ^ ^ Sin Sco s t+ °( m )• 


13 o . 68 


4 


3 4 


(3.33) 


(3.34) 


(3.35) 


Derivation of direct solar effect 

The geometry for the direct solar perturbation in the vicinity of a eollinear point is shown in 
Figure 3.3. The contribution of this direct solar effect to the effective potential energy (for motion relative 
to the Earth) is 


U D = -Afc 


r s • r 


l r s - r l 


where 


r = [(! ± XlX 1 + P) + x ]* + yi> 
r s = r s [cos (o s ti - sin co s £jJ, 

and Af@ is the mass of the Sun. Since r s » r, Equation 3.36 can be expanded in a power series 



% 

°D=~r< 1 + 

s 


3 / s 


r„ • r> 


8 V«S/ s v 


l/r_ 

SJ 




(3.36) 

(3.37) 

(3.38) 


(3.39) 


Figure 3.3-Geometry for direct solar perturbation near 
a eollinear point of the Earth-Moon system. 


Therefore, the acceleration can be written 
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P D = - V£/ D = m 2r j [ 3 (f • f s tf s - f] + (£)[f(f • *,)*, - 3(t • ^ -|f J + . . . J, 


where Af^/r 3 = m 2 . The components of P D are approximately 


P xD ~ P D * 1 = ffl2f U( 3 cos 2cu s £ + 1) ~j(y ) sin 2w s £ 


and 


+ cos (Ogt + 5 cos 3ca s t) -^^(sin co s t + 5 sin 3<y s t)j| 


P„ r, = P, 


• j = m 2 r | — ~ sin 2co J + 7 r^~)(l - 3 cos 2co„t) 


© "3 3 /y 

-(sin co s t + 5 sin 3co s £) )(cos - 5 cos Sco^t) 


>]| 


(3.40) 


(3.41) 


(3.42) 


Cancellation of first-order terms at a collinear point 


The indirect solar perturbation enters through p and v. From Equations 2.43a and 2.43b, it 
can be seen that the indirect acceleration at a collinear point is approximately 


P XI = -(l±y L Xp-2v-Sp) (3.43) 

and 

P yI -~a±y L )(2p+ul (3.44) 

To be in harmony with Colombo's solution, it is only necessary to show that all acceleration terms to 0(m 3 ) 
vanish. The linear approximations in Equations 3.43 and 3.44 are sufficiently accurate for this purpose. To 
the same order of approximation, Equations 3.41 and 3.42 are 


P xD ± y L X 3 cos 2c o s t + 1) (3.45) 

and 

P yD a -|m 2 (l ± y L ) sin 2<o s t. (3.46) 

Substituting Equations 3.31 to 3.35 [keeping only terms to 0(m 3 )] into Equations 3.43 and 3.44, it is found 
that P xI + P xD = 0 and P yI + P yD = 0 to 0(m 3 ).* 


Periodic orbit correction 

Using Equations 2.43a, 2.43b, 3.41, and 3.42 and neglecting terms of 0(m 4 ), the linearized 
equations of motion in the vicinity of a collinear point are 


*In an earlier paper (Reference 7), this author reached an erroneous conclusion regarding this cancellation because 
of inappropriate expressions for p and v. 
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As before, the first approximation to a periodic orbit around a collinear point is given by Equation 3.8. A 
frequency correction to this periodic orbit can be obtained from Equation 3.47. The constant coefficients of 
0(m 2 ) in Equation 3.47 lead to a modified characteristic equation (cf. Equation 2.49) 

s 4 - [(B l -2)+~(B l + 2)] s 2 - (2 B l + 1 )(B l - l)(l = 0. (3.48) 

For the L g point, the modified characteristic equation yields a corrected frequency of co n = 1.86247. How- 
ever, because = 1.86265, this correction is usually negligible. 

The second-order correction to the periodic orbit is found by employing the methods of Section 
B.3 of this chapter. Although the expressions for p and v in Equations 3.31 and 3.32 contain periodic terms 
with frequencies co s and 2w s when quantities of 0(/n 3 ) are retained, only terms of frequency 2 co s are consid- 
ered in this section. This should be sufficient for an order of magnitude estimate of the correction. The 
differential equations for the solar perturbation correction are 

- ( 2B L + D x „2 = K x sin + K 'x sin K - 

and (3.49) 

y n2 + 2 *n2 + ( B L - V>y n2 ~ K y cos ( 2co s + + K y cos ( w ji ~ 2o> g )£, 

where 

K k = jk[vg - 3 B lPs + f m 2 ] - +i[v 8 -|» 2 ] Ji» 7l . (3.50) 

K'x = jk[v, - 3 B lPs +| ffl 2 ] - v s oi n -|[*„ -f» 8 ] [V (3-5D 

K V + 1” 2 ) ‘ V 3 a n\ + ["3 +Y B lP 3 (3 ' 52) 

and 

K ‘y = j- {| ("3 +1“ 2 ) + -3%] + [>'3 + 1 S LP3 “I pyl' (3-53) 

with 

Ps s - (“ a+J f' nl3 )’ 
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and 


_ /II 2 , 85 3 ] 
V 3-( 4 m + 12 m 


/II 2^26 3, 


Vo = 


and 
Once again, 


\ / 

Particular solutions for these differential equations are 

X n2 = A X2 sin + w n ) £ + A 'x2 sin 
y jj 2 ~ A y2 (2co g "t C0jj)t "t -^yg cos (<y^ — 2c0g)t. 

in, consider the example of Section B.3 of this chapter. Using the values 

p 3 = -0.00692059, = 0.0183515, 


-0.0344011, 

0.07480133, 


(O c 


it is found that 


0.925199, 

B L2 = 3.19042, 
co n = 1.86265, k = 0.343336, 

x = -0.025096 A yV A x2 = 0.0029012 A v 

x = 0.017698 A yV A' x2 = -0.0024375 A yV 


K y = -0.035165 A yV A y2 = 0.0048903 A yl 


A y l ’ 

A Vl’ 


k; 


and for A yl = 0.02 = 7688.1 km 


y = -0.026236 A yl , A yQ = -0.011951 A v 


K x = ~ 1.3934 x 10" 7 g, A x2 


22.30 km, 


K' x = 9.8264 x 10 8 g, A' x2 = -18.74 km, 


K y = -1.9524 x 10"' 1 g, 
K y = -1.4567 x 10" 7 g, 


>2 


37.60 km, 


/ - y2 

Although these corrections are small, the solar perturbation terms 
culations because accelerations of 10" 7 g are not insignificant. 


(3.54 


4,o = -91.88 km. 

should be included in 
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2. Effect of the Moon at the Sun-Earth Collinear Points 


The geometry for the Sun-Earth-Moon system is depicted in Figure 3.4. Consider a point 0 that is 
collinear with the Earth-Moon barycenter B and the Sun. The distance r BQ is determined by assuming that 
the total mass of the Earth and Moon are concentrated at B, and then taking r BQ as the libration-point dis- 
tance y L for fi ~ (m 1 + + m 2 + m g ). Does the point 0 lie closer to an equilibrium solution in the 

Sun-Earth-Moon system than the classical Sun-Earth collinear libration point? This question can be re- 
solved by finding the perturbed path (dynamic equilibrium point or “nominal path”) for a satellite in the 
vicinity of the point 0. If the maximum radius of the perturbed path relative to 0 is smaller than the dis- 
tance from 0 to the classical collinear point, then it is more appropriate to consider 0 to be a collinear 
libration point of the Sun-Earth-Moon system. 

The accelerations relative to an inertial reference system can be obtained with the aid of Figure 
3.4. Defining r i; - = r ; - r if these accelerations are 


ra i r Oi , m 2 r 02 , ra 3 r 03 
f 0 - "73 + “T3 + ,3 


01 


02 


03 


(3.55) 


m 2 r l2 ffl 3 r 13 
r l ” ,3 + ,3 


12 


13 


(3.56) 


and 




m l r 12 

, ffl 3 r 23 


' r 2 

” r 3 

f 12 

r 3 ’ 
f 23 



m l r l3 

m 2 r 23 


r 3 

3 

r 13 

CO 
CO <M 

i 

of the barycenter (m 

l r Bl 

+ m 2 r B2 = 

0), it is 

r B = r 2~ l B2 ” r 2 

m ^ 

m l 

+ m 2 (r 2 

r i) “ 

1 m j 


m < 


m. 


r n +■ 


(3.57) 


(3.58) 


(3.59) 


m, (EARTH) 



Differentiating Equation 3.59 and substituting Equa- 
tions 3.56 and 3.57 gives 


m< 




1*13 


m 2 r 23 


( ffl l + ra 2) r 13 ( m l + m 2>23 


(3.60) 


Figure 3.4-Geometry for a collinear point in the Sun- This relation, along with Equation 3.58, leads 
Earth-Moon system. to* 

*This result is obtained by keeping only the first term in an expansion where the ratio of the second term to the 
first term is about 8 x 10” 8 . 
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(3.61) 


m l + m 2 + m S 


r B 3 r 3 *B m i + m < 


m l r l3 | m 2 r 23 


13 


23 


(m 1 + m 2 + ra 3 ) 


-r 


S3 1 


S3 


Making use of Equations 3.55, 3.60, and 3.61, it is now possible to give the acceleration of a satellite at 0 
relative to the Sun as 


r 30 ~ r FO r F3 ~ r O r s r B3 


m 3 r 30 , ( m i + m 2 + m g) 1/ 

~F~ + 3 r B3 ~JT (“‘Z'BZ + “l r B 0 ) 

x 30 r B3 01 


m , 


m, 


+ i3 ( r B2 r BO ) 


,3 ( m l r B3 + m 2 l B2^ 


m s m 2 


r 3 ^ r B3 r B2^- 


(3.62) 


1 02 (7^1 + /Hg)f^g ^ ^ + 

Neglecting terms of 0(r B2 /r BO ) 4 , (r B2 /r B0 = 1/4), a moderate amount of algebra reduces Equation 3.62 to 


where 


30 


ffl 3 r 3 0 


30 


(^1 + m 2 ) 


l B0 l B3 


BO 'S3 


+ P, 


m e 


BO 


r B 2 ' r B0 3 ( r B2 


r 2 

r B0 


2 \r 


BO) 


15 j B2 ‘ r B0 ] 

2 l 4o J 


r B2 + 


3 / r B2 


2\r 


BOj 


(3.63) 


15 | B2 ' r B0 | 15 / r B2] / r B2 ’ t B() \ _ 35 / r B2 ‘ T Bo f 

2 \ r§. ^ / 2 \r on) V r? _ / 2 1 »-2 


BO 


BO/ 


BO 


BO 


B07 


(3.64) 


Since r fi Q = y^, it is immediately obvious from Equation 3.63 that r 30 = P. 

For the remainder of the analysis, the orbital eccentricities and the inclination of the Earth-Moon 
orbit to the ecliptic are neglected. Essentially nothing is lost here because this approximation is con- 
sistent with the deletion of terms of 0(r B2 /r BO ) 4 . Using t B2 • r BQ = +r B2 r B 0 cos 0 (upper sign for L 2 , 
lower sign for L/), the linearized equations of motion relative to 0 become 


and 


where 


x - 2 y - (2 B l + l)x = P x 
y + 2x + (B L - 1 )y = P y , 


/^BP \ 1 2 1 

+-r[l + 3 cos 26] +ir[3 cos (9 + 5 cos 30] 

r B0\B0l ‘ 4 2 


(+J 


(3.65) 


(3.66) 


3 2 r B 2 


1 


B2 


P i — - | /+ sin 26 - ^[sin 6 + 5 sin 3(9] ' 

7 2 r BOX BO/ 4 \ r BOj 


(3.67) 
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and 6 - nt (n - 13.369). A particular solution of Equation 3.65 is the perturbed path 


and 


x n = x c + A x cos 6 + A' x cos 26 + A' x cos 3 6 


y = A v sin 6 + A' sin 26 + A" sin 30. 


(3.68) 


For the perturbed path at the L 2 point, 


m 2 = 3.7397 x 10 


r BQ = 1.0078 x 10 2 


b L2 = 3 - 9405 > 


= 0.25187, 



K = = 2.5358 x 10 5 = 1.4120 x 10 8 g, 

r 2 U 


x c = 295.10 km, 


A x = -7.47 km, 


A y ~ 3.01 km, 


^ = 11.47 km. 


A 1 = -8.22 km, 


A' x = -1.41 km. 


Ay - 1.10 km. 


The constant term x c represents a shift along the x-axis to a new point O . Relative to O , the maximum 
distance to the perturbed path is about 20 km. Therefore, the point O' is never further than about 20 km from 
a dynamic “equilibrium” point of the restricted four-body problem (Sun, Earth, Moon, and satellite). 

3. Jupiter’s Effect at the Sun-Earth Isosceles-Triangle Points 

A satellite moving in the Earth’s orbit around the Sun will be perturbed by several bodies, e.g., 
Earth, Moon, or Jupiter. The relative importance of the different effects may be quite variable, depending 
on the distance of the satellite from the Earth. In Chapter II, the Earth’s effect was considered, and artifi- 
cial libration points were created by canceling this effect with a constant thrust acceleration. The principal 
objective of this section is to obtain an estimate of Jupiter’s effect at these artificial libration points 
(isosceles-triangle points) when the satellite is relatively close to the Earth. 

The geometry for Jupiter’s perturbation at a Sun-Earth isosceles-triangle point is shown in Figure 
3.5. Neglecting some small indirect effects, Jupiter’s perturbative acceleration is 




■_E_ 

r 3 

EJ 


(3.69: 


Assuming that (a/r) is small, Equation 3.69 becomes 
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Figure 3.5— Geometry for Jupiter's perturbation at a 
Sun-Earth isosceles-triangle point. 



and the components are given by 

p jx = p j- 1= l5(r) sin2a (3 - 71) 

and 

p iy = p j ■ j = i p(t ) [1 + 3 cos 2al (3,72) 

The magnitude of Jupiter's acceleration is 

i p ji =| 7(r) [1 ° + 6 cos 2a]1/z - (3 - 73) 

For 0 = 10° and a - 90°, rough calculations give 
|p | =1.35 x l0~ 9 g (this value is relatively constant 
for 0° < a < 180°). Because the thrust acceleration 
at 6 =10° is |F C | = 6.02 x 10~ S g, Jupiter’s pertur- 
bation can be ignored when 6 < 10°. 


4. Solar Effect at the Earth-Moon Equilateral-Triangle Points 

The analysis of the uncontrolled motion of a satellite in the vicinity of the Sun- perturbed Earth- 
Moon equilateral-triangle points represents one of the most difficult problems in celestial mechanics. In the 
past few years, several papers, e.g., see References 48 to 55, have appeared on this subject, A recent ana- 
lytical study by Schechter (Reference 54) is probably the most productive to date, and some of his conclu- 
sions have been verified by the numerical work of Kolenkiewicz and Carpenter (Reference 55). However, 
the general problem remains unsolved because the important effect of the lunar eccentricity was not consid- 

ered in Schechter 's study. 
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CHAPTER IV 


LINEAR FEEDBACK CONTROL 


After a satellite has been injected into the vicinity of a libration point, a translation control system can 
be used to eliminate the initial errors. This control system can also be used for station keeping. To imple- 
ment this control, closed- loop guidance laws (feedback control logic) will be needed. In this chapter, sim- 
ple linear feedback control laws are developed for the collinear and equilateral-triangle points. Although 
the control acceleration is treated as if it were continuous, this assumption is not overly restrictive because 
a pulsed control will behave essentially like a continuous control if the frequency of the pulse control is 
much higher than any of the natural frequencies of the system.* 

For motion in the immediate vicinity of a libration point, the linearized equations of motion are suffi 
ciently accurate. These equations will contain periodic coefficients if eccentricity is present or if the 
satellite is following a nominal path. In many instances, these effects are small and can be neglected. 
However, to obtain some estimate of the importance of these periodic coefficients, Floquet stability investi- 
gations for some representative cases are included here. 

A. Collinear Points 

For small eccentricities, the linearized equations of motion at a collinear point (or relative to a small 
nominal path about the point) are approximately (from Equation 3.4)t 


x - 2y - (2 B l + l)x = F cx , 

(4.1a) 

y + 2x + (B l - l)y = F cy , 

(4.1b) 

Z+ B L z = F cz- 

(4.1c) 


The z-axis motion is bounded and can be damped by simply taking F cz = - k^z, where the feedback gain 
k ' > 0. The control synthesis for the unstable xy - motion is less obvious, but it can be accomplished with 
classical methods [Routh criterion (Reference 56) and root locus (Reference 57)]. 


*The continuous approximation is probably adequate if the pulsing rate is greater than 10 times any of the natural 
frequencies of the system. 

tFor the equations relative to a nominal path, the variables (g, p, Q should be used {g-x- x n , etc.). However, 
this distinction will be overlooked when the equations have the same form. 
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1. Routh Stability Conditions 

For a radial-axis (x-axis) control using only range and range-rate (x and x) feedback, the control is 


and 


F cx = -*1* - V 


F.y = 0 . 


(4.2) 


This leads to the characteristic equation 

s 4 + k t s 3 + [k 2 - ( B l - 2)]s 2 + k 1 (B L - 1 )s + [kg - (2 B L + 1 )](S L - 1) = 0. 
From Routh’ s criterion, the necessary and sufficient conditions for asymptotic stability are 


(4.3) 


k 1>0 


and 


kg > (2 B l + 1). 


(4.4) 


It is shown in the next section that, for certain values of k 1 and kg, more damping can be obtained by adding 
some positive y-feedback. For this case 


F cx = ~ k l *- v + v 


and 


F cy = 0. 


(4.5) 


and the characteristic equation is* 

s 4 + k 1 s s + [ig - (B l - 2)]s 2 + [tj(B L - 1) + 2* 4 ]s + [*g - (2B l + 1)](S L - 1) = 0. 


(4.6) 


The Routh conditions are 


and 


kj >0, kg > 2B l + 1, 
k x (B l - 1) + 2k 4 > 0, 

2k~(B L - 1) - Jc 4 |t 1 [3B L - (k 2 + 4)] + 2k 4 |> 0. 


(4.7) 


If cross-axis (y-axis) control is used, it is necessary to feed back all of the state variables for 
stability. This control can be written 


F = 0 

cx 


and 


F cy ~ - V " V " V + V- 


(4.8) 


*The same characteristic equation is obtained with the dual-axis control, F cx = ~ k 1 x - k 2 x, F cy = -k 6 x, with k e 
replacing k 4 . 
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The characteristic equation is given by* 


s 4 + k 7 s s + [ 2 fr 5 ~ k s - (. B l - 2)]s 2 + [2k 6 - k 7 (2B L + 1 )]s + [k g - (B L - 1 )](2B L + 1) = 0, (4.9) 

and the Routh conditions are 

k 7 > 0, 2k & > k s + (B l - 2), 

2k & > k 7 (2B L + 1), k 8 > (B l - 1), (4.10) 

and 

fc 7 [2fr 5 ~ k 8 + (3 B l + 4)] - 2k A ~ k*(k 5 + 2)(2B L + 1) > 0. 

2. Root Loci and Closed-Loop Response 

The preceding section has shown that it is possible to obtain a stable radial-axis control that re- 
quires only range and range-rate measurements. However, it is still necessary to select the gains k 1 and k 2 
for a satisfactory response. A simple method for determining these gains is to use root-locus plots. Taking 
k 1 = 0 and writing Equation 4.3 in the form 

fc a [s g + (B L -l)] 

s 4 - (B l - 2)s 2 - (2 B l + 1)(B l - 1) “ ” 1 (4 ' U) 

leads to the root-locus plot of Figure 4.1. It is seen 
that neutral stability is obtained when k 2 >(2B L + 1). 
When k 1 ^ 0, Equation 4.3 can be written 

+ (B l - 1)] 

s 4 + [k 2 - (B l - 2)]s 8 + (B l - 1)|V (2B l + 1)] (4 12) 

The root loci of Equation 4.12 and some correspond- 
ing response diagrams! are depicted in Figures 4.2 
to 4.7 for three values of k 2 . Inspection of these 
figures shows that, for B L = 4, a well-damped re- 
sponse can be achieved by taking k^ - 4.0 and = 
12.0 (note that k 1 =- 4 -k g ). 

It is obvious that the response will be uh- 
derdamped whenever k 2 is too large. However, when 
the radial-axis control is augmented with positive 
y-feedback, it is possible to obtain more damping. 

*This characteristic equation can also be obtained with the dual-axis control, F cx = k%y + Jr 4 y, F cy = -k 7 y+ k 8 y, 
with lc 3 and k 4 replacing k 5 and k 6 , respectively. 

|The response diagrams were obtained from a TR-48 analog computer simulation. 
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Figure 4.2-Left, root-locus plot for 
Equation 4.12. Gain is ky k 2 = 10, 
and B L = 4. 


Figure 4.3-Below, response with 
radial-axis control of Equation 4.2, 
y versus x, k 2 ~ 10, B L - 4, and 
x(0) = y(0) = 0 . 
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Figure 4.4-Left, root-locus plot for 
Equation 4.12. Gain is k ^ k 2 - 12, 
and B l - 4. 


Figure 4.5-Below, response with 
radial-axis control of Equation 4.2, 
y versus x, k 2 = 12, B L = 4, and 
*(0) = y(0) = 0. 
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Figure 4.6-Left, root-locus plot for 
Equation 4.12. Gain is k p k 2 = 
16, and B L = 4. 


Figure 4.7-Below, response with 
radial-axis control of Equation 4.2, 
y versus x, k% - 16, B L - 4, and 
x(0) = y(0) = 0. 
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This can be seen at once by comparing Figures 4.7 and 4.8 (note instability at k 1 - 0.5 in Figure 4.8.) An 
explanation for this behavior can be found by writing Equation 4.9 as 

2k * S F _=-l (4.13) 

s 4 + kjS 3 + [kg - (B l - 2)]s s + k j(B l - l)s + (B L - l)[k 2 - (2B L + 1)J 

and examining the resulting root-locus plot given in Figure 4.9. It can be seen that the damping will be 
Increased for small values of k 4 , but instability can occur if k 4 becomes too large. 

3. Floquet Stability Investigations 

In some instances, the magnitudes of the periodic coefficients in the linearized equations of motion 
are appreciable. Therefore, the effects of these coefficients on the stability of the closed-loop system 
should be examined. If the effects from the eccentricity and the amplitude of a nominal path are analyzed 



k 1 = 4.5. 

Figure 4.8-Response with radial-axis control of Equation 4.5, y versus x, 
k 2 - 16, = 3.5, B l = 4, and x(0) = y(0) = 0. 


53 



Figure 4. 9-Root- locus plot for Equation 4.13. Gain 
is k 4 , k l = 4.5, k 2 = 16, and B L = 4. 


separately,* this examination can be performed in an 
efficient manner by making use of Floquet theory 
(Reference 58). 

The linearized equations of motion can be 
written in the matrix form 

X(t) = F(t)X(t). (4.14) 

For a two-dimensional analysis (fourth-order sys- 
tem), X(t) is a 4 x 1 state vector and F(t) is a 4 x 4 
periodic matrix [i.e., F(t) - F(t + T)]. To implement 
Floquet’ s theory, it is first necessary to define a 
4x4 matrix H(t) by 

H{t) = F{t)H(t) 

and (4.15) 

H( 0) = / (unit matrix). 

Equation 4.15 is integrated numerically to obtain 
H(T), and the following characteristic equation is 
formed 

det [ H(T ) - A/] = 0. (4.16) 


On the assumption that the roots A ■ (j = 1, 2, 3, 4) of Equation 4.16 are distinct, Floquet’s theorem states 
that all solutions of Equation 4.14 are bounded if and only if 

|Aj.|<1 j = 1,2, 3,4. (4.17) 


Eccentricity effect 

From Equations 2.43 and 4.2, the linearized equations of motion at a collinear point with 
radial-axis control are 

x - 2(1 + u)y - ( 2B l + l)x = (2v - 6B L p)x + i>y - k^x - k g x 

and (4.18) 

y + 2(1 + v)x + (Bj^ - l)y = - vx + (2v + 3B L p)y. 

Neglecting terms of 0(e 2 ) (cf. Equation 3.19), 

p - -e cos t 

and (4.19) 

v - 2e cos t. 

Therefore, T = 277, and 


*In general, some coupling is present. 
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and 


*(*) = 


x(t) 

x(t) 

y(t) 

m 


F(t) = 


0 

v 

' 21 

0 

^41 


1 

-k 

0 


42 


0 

23 

0 

^43 


24 1 

1 

0 


F 21 “ \( 2B L + !) - + 2e ( 2 + 3J3 l) C0S *]’ 


F 23 = ~ F 4i - sin t, . 


F 24 = “ F 42 = + 2e C0S ^ 

F 43 = - [( B L - 1) + e (3B L - 4) 008 l ] 


(4.20) 


(4.21) 


(4.22) 

(4.23) 

(4.24) 

(4.25) 


When Equations 4.20 to 4.25 are substituted into Equation 4.14, and when the aforementioned procedure is 
carried out, the instability regions of the parameter space (k^, k 2 ) are determined. Instability charts for the 
collinear points of the Earth-Moon system are given in Figure 4.10. 

Two interesting features of these charts should be pointed out. By analogy with the constant 
coefficient system, Equation 4.22 suggests that the periodic coefficient system may be stable if 

k 2 > (2 B l + 1) + 2e(2 + 3 B L ). (4.26) 

For the Earth-Moon system, Equation 4.26 gives > 13.19 and (k > 8.65. When these inequalities 

are satisfied, the larger instability regions of Figure 4.10 are excluded. The remaining points of instability 
can be explained by resonance. The value of k 2 for closed-loop resonance of the constant coefficient sys- 
tem can be found from Equation 4.3. For k^ = 0, this value of k 2 is given by 


(2 B l + 1 )(B l - 1) - (B l - 2 )co 2 - v 4 

(B l - 1) - CO 2 


(4.27) 


When a) = 1, the critical values of k 2 in the Earth-Moon system are (fc 2 ) L1 - 13.57 and (k<^ L2 - 11.74. Figure 
4.10 shows that points of instability are located in the vicinity of these critical values. 


It is tempting to assume that Equations 4.26 and 4.27 can be used to find the principal insta- 
bility regions for any system. At the L 2 point of the Sun-Mercury system, Equation 4.26 gives k 2 > 14.68, 
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Figure 4. 10- Instability charts for controlled motion, with radial- axis control of Equation 4.2, at the 
collinear points of the Earth-Moon system, e = 0.05490, and © - unstable. 


and Equation 4.27 yields k 2 = 11.98. A comparison of these values with the instability chart of Figure 4.11 
shows that a more sophisticated interpretation is needed for larger eccentricities. 


Motion relative to a periodic orbit 


If eccentricity terms are neglected, the linearized equations of motion relative to a periodic 
orbit around a collinear point with radial-axis control are (see Equations 3.30 and 4.2) 


and 


’£ " 2 t) - ( 2 B l + 1 )£ = + 3C L 4 yl [p(sin <u n t)£ - (cos <u n t) 77 ]- kj - k^ 

'ij + 2£+(B L - 1)t/ = ± 3C L 4 yl |k(sin co n t)rj + (cos &> n t)^] . 

In this case, T = 2n/co n , and 


X(t) 


m 

At) 
v (0 
i(t) 


F(t) 


0 

7 

21 

0 

7 

41 


1 

‘*1 

0 

-2 


0 

1 ? 

23 

0 

'43 


F 21 = 


[(2£ l + 1) - 


k 2 + 6C L kA yl sin <u 




(4.28) 


(4.29) 


(4.30) 


(4.31) 
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^23 ~ ^41 ~ 1 C0S 


(4.32) 


Figure 4. 11-Instability chart for controlled motion, 
with radial-axis control of Equation 4.2, at the Lg 
point of the Sun-Mercury system. B L2 - 3.9772, e - 
0.20563, and© - unstable. 


and 


F 43 = - 1(B l - !) + 3C L A-4 yl sin <n n t], (4.33) 


The instability regions of the parameter space (A yl , 
k^, k 2 ) can be obtained by substituting Equations 
4.29 to 4.33 into Equation 4.14 and following the 
procedure given earlier. Instability charts for the 
L g point of the Earth-Moon system are given in Fig- 
ure 4.12 for six values of the gain For a lightly 
damped system (k 1 = 0.1), the region of instability is 
rather large. However, the size of this region 
shrinks quite rapidly as the damping is increased. 

At kj = 2.5, it has almost disappeared from the param- 


eter space of Figure 4.12. For co = <*> n /2 = 0.931325, 
Equation 4.27 gives kg = 10.87. Thus, the insta- 
bility region of Figure 4.12 is caused by a subharmonic resonance at a> = &> n /2 as might be expected from 
the form of F(t). 


B. Equilateral-Triangle Points 

If eccentricity is neglected, the linearized equations of motion at an equilateral-triangle point are 
(from Equation 2.62) 


X' - 2y' - ax' ~ F cx t, 

(4.34a) 

y' + 2x' - (3y‘ = F cyl , 

(4.34b) 

z ' + z ' = F cz" 

(4.34c) 


Once again, the z'-axis motion can be damped by using the control F cz , = - k^z' with k^ > 0. As shown in 
Chapter II, for p < 0.03852 . . . , the coupled motion in the x'y'-plane is neutrally stable when F cx , -F , =0. 
In this section, a simple single-axis feedback control is designed; this control will provide adequate damp- 
ing for the oscillatory motion.* 

1. Routh Stability Conditions 

Damping can be accomplished with the radial-axis (x'-axis) control 

F cx , = - fc,x' + k 4 y> 

and (4.35) 

F cy , = 0. 


♦Fleming (Reference 5) has designed a dual-axis feedback control to perform the same task. 
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*1 = 1.5. 




*! = 2.0. *1 = 2.5. 

Figure 4. 12-Instability charts for controlled motion, with radial- axis control of Equation 4.2, relative to 
a periodic orbit around the L 2 point of the Earth-Moon system. © —unstable. 


With this control, the characteristic equation becomes* 

s 4 + fcjS 3 + s 2 + (2 k A - k^)s + a/3 = 0, (4.36) 

*The dual-axis control, F cx , = ~k t x', F cyt = - k 6 x ' , yields the same characteristic equation, with k 6 replacing k 4 . 
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and the Routh conditions are 


> 0, 2k A > & l (3, 

and (4.37) 

*![* 4 (1 + 2/3) - 2ft i/3] > 2it|. 

The damping can also be obtained with the cross-axis (y'^axis) control 

F cx' = 0 

and (4.38) 

F cy' = - k 6 x ' - V- 

For this case, the characteristic equation is* 

s 4 + k ? s 3 + s 2 + (2k Q - k 7 a)s + a/3 - 0, (4.39) 


and the Routh conditions are 


k 7 > 0, 2 k Q > k 7 a, 

and (4.40) 

k^k G (l + 2a) - 2k ? a]> 2 Jr|. 

Notice that Equation 4.40 has the same form as Equation 4.37, with a replacing (3. 

2. Root Loci and Closed-Loop Response 

Because of the similarity of Equations 4.36 and 4.39, the root-locus plots for both the radial-axis 
and the cross-axis controls are equivalent. Therefore, only the radial- axis control of Equation 4.35 is con- 
sidered here. To obtain quantitative results, the parameters a and (3 (see Equations 2.63 and 2.64) are com- 
puted for the Earth-Moon system, and 


a = 2.97275, /3 = 0.02725. 


(4.41) 


If k 4 = 0, Equation 4.36 can be written as 

k^s 2 - /3) 
s 4 + s 2 + a/3 


(4.42) 


and the root-locus plot of Figure 4.13 is obtained. By choosing appropriate values for k ^ and k 4 , the unsta- 
ble roots are taken into the left half of the plane. This can be seen graphically by writing Equation 4.36 as 


s 4 + k^s s + s 2 - k^fis + a/3 


(4.43) 


*The dual-axis control, F cx i ~ k 4 y ' , F cy i = - k 7 y ' , yields the same characteristic equation, with k 4 replacing k 6 . 
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and inspecting a root-locus plot for a fixed value of k v For k t = 2.0, the root-locus plot of Equation 4.43 
is given in Figure 4.14. The Routh conditions of Equation 4.37 show that asymptotic stability is obtained 
when 

k A 

0.05808 <~r~ < 0.4692. (4.44) 

In Figure 4.15, the closed-loop response for k 1 = 2.0 and k 4 = 0.25 is compared with the uncontrolled motion. 



Figure 4.13-Root-locus plot for Equation 4.42. Gain Figure 4.14~Root-locus plot for Equation 4.43. Gain 
is k v a = 2.97275, and £ = 0.02725. is k 4 , k t = 2.0, a = 2.97275, and /3 = 0.02725. 
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Figure 4.15-Response at the equilateral-triangle points, y' versus x'. a - 2.97275, /3 0.02725 

and x'(0) = x'(0) = y '(0 ) = 0. 


Without control . 


With radial-axis control of Equation 4.35, 
k l = 2.0 and lr 4 = 0.25. 



CHAPTER V 


STATION KEEPING 


Several methods for satellite station keeping in the vicinity of a collinear libration point are examined 
in Chapters V to VII. An on-off control system is analyzed in Chapter VI, and an unconventional method is 
investigated in Chapter VII. In the present chapter, the radial-axis control of Equation 4.2 is treated. Gen- 
eral analytical relationships for the control requirements are formulated, and a solar sail control technique 

is also presented. 

As noted earlier, station-keeping problems for libration-point satellites have not received very much 
attention. On the other hand, station-keeping techniques for synchronous satellites of the Earth (24-hour 
satellites) have been thoroughly analyzed (References 59 to 63), and a large reservoir of practical experi- 
ence has been acquired (References 64 to 66). This knowledge could be very useful in the design of propul- 
sion systems, sensors, attitude control systems, etc., for libration-point satellites. 

A. Average Control-Acceleration Requirements 

Realistic estimates of station-keeping costs are needed for feasibility studies of libration-point satel- 
lite missions. In this section, some general expressions for these costs are derived. To establish a basis 
for comparison, note that the average control accelerations for synchronous satellites are about 1.5 x 10 g 
(1 fps/yr = 1 x 10~ 9 g) for north-south station keeping, and about one-tenth of this value for east-west station 

keeping (Reference 63). 

1. Cost Estimates for Noise Inputs 

The station-keeping cost for a satellite that is following a nominal path around a collinear libration 
point is a function of measurement noise, engine fluctuations, and random accelerations from such causes as 
gas leakage or random solar radiation pressure effects. A block diagram of the radial-axis control system 
with noise inputs is shown in Figure 5.1. The measurement noise is denoted by n* n x ; the engine noise is 
n ex ; and the random accelerations are n px , n py . A simple lag filter has been added to reduce the effect of 
the measurement noise. An optimum design can be achieved by choosing the control parameters ( k v k 2 , r) 
to minimize the mean square value of the control acceleration F cx (a description of this minimization proce- 
dure can be found in Reference 67). Because the uncontrolled system is unstable, this minimization is 
meaningful even when n px = n py = 0. 
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Optimization 


P* 



From Figure 5.1, it can be seen that 
the control acceleration is 


F cx(s) TS 


+ k 2 )x ^ 
+ (tji + k s n x )J + n e 


(5.1) 


which can be written as 
(Vx + k 2 n x> 


=■ 


h(s) 

(2B l + 1 )(B l - 1) 
(■ k ! s + k 2 ) 


[® 4 - ( b l - 3)s 2 

>] 


Figure 5.1-Radial-axis control system with noise 
inputs. 


h(s ) 


) 


s 8 + (B, - 1) I (n - n ov ) 


+ 2s V + n ax, 


px ex' 

(5.2) 


where 


A(s) = rs 5 + s 4 +^k x - t{B l - 2)Js 3 +^ 2 - (B l - 2)Js 2 

+ ^ - r(2B L + l)](B L - l)s + [kg - (28 l + 1)](Z? L - 1). 


(5.3) 


For the characteristic equation, h(s ) = 0, the Routh stability conditions are 

r > 0, 

k 2 > (2 B l + 1), 
and 

Writing the characteristic equation in the form 

/ k r 


(5.4) 


s + T 


S 2 + (B l - 1) 


(rs + l)[s 4 - (B L - a )s 2 - (2 b l + 1 )(B L - 1)] 


= - 1, 


(5.5) 


the root-locus plot of Figure 5.2 is obtained. 

The engine noise and random accelerations can usually be neglected in a first approximation, 
and the control acceleration then becomes 


64 








(5.6) 


F cx( s ) = 


^l n x 4 k 2 n x ) 


h(s) 


s 4 - (B l ~ 2 )s 2 - (2 B l + 1 )(B l - 1) 


Assuming that n ^ and n x are independent, white, Gaussian processes, the mean square value of Equation 
5.6 is 


r+jca 

( F cx) 2 = N i [ k l + *2 K J»] 2^ / 

~'—jco 


g(S ) 


h(s)h(~s) 


ds,iV.P 5 , 


(5.7) 


where 


g(s) = s 8 - 2 (B l - 2)s 6 + [(e L - 2) 2 - Z(2B l + 1 )(B L - l)]s 4 
+ 2(B l - 2)(2 B l + 1 )(B l - l)s 2 + (2 B l + 1 ) 2 (B l - l) 2 . 


(5.8) 


and N x are the noise levels, and R N = (iV x /l\y. The integral of Equation 5.7 can be written as an alge- 
braic function of the control parameters Ifg, r) by using a special table of integrals (Reference 67). 
Because the resulting expressions are rather cumbersome * a numerical search procedure (Reference 68) is 
used to find a minimum value of the cost parameter P 5 , while the stability conditions of Equation 5.4 are 
satisfied. It is found that the optimum values are 


Q, 


k 2 = ^ 2B L + 




(5.9) 



where Q is a function of B L and is plotted in Figure 
5.3. Therefore, the minimum average control accel- 
eration is given by* 



(5.10) 


The effects of periodic coefficients 
(see Chapter IV) have not been considered in the 
foregoing optimization. This omission is probably 
not too serious as long as the unstable regions of 
the parameter space are avoided, but some method of 
assessing the influence of these coefficients would 
be useful. Equation 4.26 implies that an effective 
value of B l , 


Figure 5.2— Root-locus plot for Equation 5.5. Gain 
is k v Not to scale. 


B' l 5E B l + e(2 + 3P l ), (5.11) 


*For a Gaussian process, F< 


( F 4 ms = 0.8(F cx ) rr 
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Figure 5.3-Cost parameter Q versus B h . 

with W = (minute)" 1 (Reference 69). Estimates of a x 
and 71. 


may be more appropriate for the general case. 
Although this conjecture seems to be reasonable, 
verification is required. 


Estimation of noise level 

The noise level can be estimated 

by using 


»i = 



(5.12) 


where a? is the range-rate variance, W is the sam- 
pie rate, is the time interval for a single set of 
measurements, and T ^ is the time interval between 
sets of measurements. For deep space instrumenta- 
tion facility (DSIF) measurements, = 0.02 m/sec 
for onboard sensors can be found in References 70 


Costs at the Earth-Moon collinear points 

The general results of Equation 5.9 can be used to find typical values for the station-keeping 
costs at the Earth -Moon collinear points. For the examples given here, periodic-coefficient stability results 
are taken into account while choosing the ratio k^/kg, but Equation 5.11 is not employed. 

With DSIF measurements* and T M = 1 day and T 1 - 5 minutes, Equation 5.12 gives N k - 1.757 
10“ 11 = 6.912 m 2 /sec. The ratio k^/kg is chosen to be 0.1, which means that R N = 0.01. Taking N x = 
p/^s T i) T M solving for o x yields (j x = 0.751 km. (This accuracy is easily attained with DSIF meas- 
urements.) Equations 5.9 and 5.10, along with Figure 5.3, give for L 1 


- 1.505, kg = 15.05, 

t — 0.033, P 5 = 204, 

and „ 

|F cx j = 4.790 x 10 5 = 1.330 x 10 

and for L 2 


and 


k x = 0.9841, k 2 = 9.841, 
r - 0.033, P 5 = 86, 

IF - i = 3.110 x 10“ 5 = 8.858 x 10“ 9 g. 


*It is assumed here that the satellite in the vicinity of L z is following a quasi-periodic orbit around the point and 
is visible from the Earth tracking station. 
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The values of |F cx j for these examples are roughly equal to the average control accelerations that are 
needed for east-west station keeping of synchronous satellites. 


2. Cost for Sinusoidal Control Acceleration with Noise 

In some instances (see Chapter VIII), the station-keeping cost is given by |F* | where 

F cx = P( t ) + ^ cx . (5.13) 

p(t) = K c cos cot, (5.14) 


and F cx is a random variable with a Gaussian distribution. Therefore, F* x has a Gaussian density function 


$ F *(x, 0 = = exp 

Op^/Zn 


[x - p(£)] 2 ) 

" z 4 ) 


(5.15) 


with a F s (f cj ) (for the radial-axis control of Figure 5.1, a F = For this ease, the average 

control acceleration is given by 


where 



(5.16) 


(5.17) 


and F{y| is the “expected value” of a random variable y. After some rather lengthy manipulations (see 
Appendix A), Equation 5.16 becomes 



(5.18) 


where J 0 (u) and / 1 (u) are modified Bessel functions and a = k\/ 4ct|. Using well-known series expansions 
(Reference 72) for I Q (u ) and I^u), two limiting cases of Equation 5.18 can be deduced: 


Case 1, for small values of u, 


2oi 


IF* I = 

I ex I 




1 + 




and Case 2, for large values of u, 


2 K r 


CX I 


1 + 


2 Kt 


(5.19) 


(5.20) 
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3. Cost for a Constant Displacement 


A satellite is sometimes stationed at a constant distance 8 from the libration point. This displace- 
ment may be intentional, or it may be the result of a bias error in the measurements. In either case, the 
station- keeping cost will be increased. If terms of 0(5^) are neglected, this cost can be obtained from 
Equation 2.43. The special cases in which 8 lies along a single coordinate axis are illustrative. If p and v 
are neglected, these costs are approximately 

for 8 X only, 
and 

for 8 y only, 


Fey = ( B L ~ -f D L*r 

and 


F cx = - ( 2B L + 1 )S X ± 3 C L Sl - 4 D h Sl 


(5.21) 


F = F = 0; 
cy cz ’ 


F cx = tf c L s r 


and for 8 Z only. 


and 


F cx=+{ C L S % 


F C y = 0 . 


Fez = B lA “f »L*r 


(5.23) 


The costs for a satellite in the vicinity of the L 2 point of the Earth-Moon system are given in Figure 5.4. 
Notice that the contribution of the nonlinear portion of Equations 5.21 to 5.23 is extremely small for this 
range of 8. 


B. Solar Sail Control at the Earth-Moon Col I inear Points 

It has been shown that the magnitude of the control acceleration needed for station keeping in the vicin- 
ity of a collinear point of the Earth-Moon system may be extremely small. Therefore, it may be feasible to 
obtain this acceleration with a small solar sail. At first, the utilization of a solar sail to perform this task 
appears to be straightforward, until it is realized that the control- acceleration vector can never have a com- 
ponent directed toward the Sun. Furthermore, the direction of the incident solar radiation is varying contin- 
uously with respect to the xyz-coordinate system. Fortunately, these difficulties can be circumvented by 
using a method that was originally devised by Colombo (Reference 3), who demonstrated that the required 
accelerations could be obtained by forcing the libration-point satellite to follow a path that is synchronized 
with the Sun’s motion. 
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Figure 5.4~Station-keeping cost for a constant dis- 
placement from the L 2 point of the Earth-Moon system. 

Moon system is depicted in Figure 5.5. For specular 
motion (cf. Equation 3.26) are 


An analysis of a particular solar sail control 
technique is presented in this section. Colombo's 
strategy is employed, and explicit relations for the 
sail variations are derived. The analysis is limited 
to the xy-plane, and the effects of eccentricity and 
solar perturbations are neglected. Solar occultation 
periods are also neglected in this analysis. Al- 
though it is possible to use the solar sail for both 
attitude and position stabilization, only position 
control is considered here. 

1. Basic Concepts 

The geometry for a plane solar sail in the 
vicinity of a collinear libration point of the Earth- 
reflection of the incident solar rays, the equations of 


and 


x - 2 y - (2 B l + l)x = - K p cos 2 ft cos (0 + ifr) - F cx 
y + 2x + (B l - l)y = K cos 2 if/ sin (0 + ifj) - F cy , 


(5.24) 


where 6 = (o s t (a> s = 0.92519867), K p and ifr are bounded control variables (K p > 0, |^| < rr/2), and from 
Equation 3.1, 


k p = 


(5.25) 


Following Colombo’s suggestion, the libration-point satellite is required to follow a path [x Q (t), y 0 (t)], which 
is determined by finding the forced response of 


and 


x 0 - 2y 0 - (2 B l + l)x Q = -K 0 cos 6 
y 0 + 2x 0 + (B l - l)y 0 = K 0 sin 0. 


(5.26) 


The connection between the control accelerations (see Figure 5.5) is given by 


with 


and 


Fc = F 0 + K 

Fo “ F C x* + F C yj, 

F 0 = K q [- cos 0 i + sin 6 j], 

k - n? + nj- 


(5.27) 

(5.28) 

(5.29) 

(5.30) 
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Figure 5.5-Plane solar sail in the vicinity of a 
eollinear libration point of the Earth-Moon system. 



r — K i 


Figure 5.6-Solar sail acceleration diagram. 

( K n ) cos 2 (/r has a common tangent with the circle 

common tangent is at if/ - 0°, and it is obvious that 


Using linear superposition, it is obvious that the 
equations of motion relative to the path [x 0 (t), y 0 (£)] 
are 

'£~2-n-(2B L + l){=F' ci 

and (5.31) 

’t}+2 £+(B l - 1 )rj = F' Crj , 

where £ = x - x Q and rj = y - y Q . 

A graphic exposition of the solar sail con- 
trol technique can be obtained by inspecting the 
acceleration diagram of Figure 5.6. In this figure, 
the magnitudes of the control accelerations needed 
for station keeping are represented by K h (measure- 
ment bias errors), K n (measurement noise), and K d 
(perturbations, etc., which are not eliminated by 
following a nominal path), 

and 

K^K b + K n + K d . (5.32) 


It is easily seen that the required acceleration can 
be provided by suitably adjusting the control varia- 
bles K p and if/. Figure 5.6 can also be used to de- 
termine [K n 1 . From Equation 5.25, it can be seen 

\ ty max 

that this specification is equivalent to finding the 
maximum size of the solar sail. The minimum re- 
quired value of IK ) can be found by observing 
that, for this special case, the curve ]F C | = 
of radius Kj. For ^ max < 30° (sin t i max = K x /k q ), the 


(* P ) m » = 

When «A max - 30°, the analysis is more involved, but it can be shown that 


(5.33) 



Mr 

4 K ° 




(5.34) 


and the common point of tangency can be determined from 
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(5.35) 


COS iff = 




1/2 


2 . 


and 


Sail Variations with Radial-Axis Control 

For a radial-axis control, F' = 0, and Equations 5.24, 5.26, and 5.31 give 

C77 

- K p cos 2 if/ cos (6 + if/) + K 0 cos 6 = F' c ^ 

K p cos 2 if/ sin (6 + if/) - K Q sin 6 = 0. 


(5.36) 


It is convenient to use the definition 


K p = K 0 + Q 


(5.37) 


where Q p represents the variable portion of K p . The sail variations iff and Q p can be expressed as functions 
of the variables F' c ^ and 6. Equations 5.36 and 5.37 yield 


and 


F‘ c e sin 6 
tan ifr - j. _ 

K 0 b 


Q, 


= Kr 


(K 0 - F' c g)/ cos 3 if/ 
K 0 - F' c g( 1 - cos 6) 



(5.38) 


(5.39) 


and if the radial-axis control of Equation 4.2 is used, 

F W = - V- 

it 

When F^ ^ Kq, SQUEtions 5.38 d^nd 5.39 l63>d to th6 3<pproKiiii3*t6 rolcttions 

if; = -T 7 - sin 6 


(5.40) 


K 


and 


o 


e P * 


(5.41) 


3. Examples 

The variation of the area of the solar sail is minimized by choosing K Q » K v (The variation of ifr 
is also minimized.) It is quite possible that the mechanization of the solar sail control system will be less 
complicated if only small sail variations are required. Therefore, a value of K Q = 15 K 1 0/f max = 3.82°) is 
adopted for the examples given below.* 

*Although sail variations are minimized when K 0 ~^> K^ t the required sail area A s is increased. This tradeoff is 
not considered here. 
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For both examples, it is assumed that K d ~ 0, 8^ - 1 km, and m - 400 kg. K n is obtained from 
section A.l of this chapter, and p B - 4.50 x 10~ 6 newtons/m 2 . The forced response of Equation 5.26 is 
written as 


x 0 “ ^xo cos co s t 


and (5.42) 

y 0 = 4 y0 sin <a s t. 

With the values mentioned above, some of the important parameters for satellites in the vicinity of the Earth- 
Moon collinear points are found to be 

for L 1 

K b = 2.938 x 10~ 5 = 8.157 x 10 ~ 9 g, 

K n = 4.790 x 10~ 5 = 1.330 x 10 _8 g, 

K Q = 1.159 x 10' 3 = 3.218 x 10 " 1 g, 

A x0 - 14.8 km, A y0 - 143.7 km, 

A s = 149.7 m 2 , and (A s /m) = 0.374 m 2 /kg, 


and for L g 

K b = 1.920 x 10“ 5 = 5.330 x 10" 9 g, 

K n = 3.110 x 10" 5 = 8.858 x 10” 9 g, 

K q = 7.545 x 10" 4 = 2.128 x 10 “ 1 g, 

* 

A x o = - 10.4 km, A yQ - 202.9 km, 

A s = 97.4 m 2 , and (A s /m) = 0.244 m 2 /kg. 

Notice that the forced oscillation is counterclockwise for L 1 and clockwise for L g . In both cases, ji4 Q 1 » 
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CHAPTER VI 


ON-OFF CONTROL SYSTEM 


It is quite possible that the station keeping of a libration-point satellite will be accomplished with an 
on-off control system. Accordingly, an analysis of the ensuing limit-cycle motion is in order. Two methods 
of analysis are presented in this chapter. The first is an exact method that makes use of Fourier expan- 
sions; the second is an approximate phase-plane technique. 

Once again, a radial-axis control that requires only range and range-rate measurements is utilized. If 
eccentricity, perturbations, etc., are neglected, the equations of motion in the vicinity of a collinear libration 
point are 


and 


where 


x- 29- ( 2B l + l)x = - F(£) 
y + 2x + ( B l - 1 )y = 0, 


£ = x + Ax, 


( 6 . 1 ) 


( 6 . 2 ) 


and the assumed on-off control characteristic F(£) has a deadzone £ d and a hysteresis loop of width 25 (see 
Figure 6.1). An asymmetric characteristic was chosen because only “one-sided” limit cycles are considered 
here.* The constant A. in Equation 6.2 is just the slope of the “switching line” in the x-x phase plane. 


A. Limit Cycles: Exact Analysis 

Conditions for the existence of limit cycles in Equation 6.1 can be found by using harmonic methods. 
Approximate methods (References 73 and 74) that neglect all harmonics higher than the first are very popu- 
lar, but their validity is sometimes questionable. Therefore, the higher harmonics are retained in the analy- 
sis employed in this section. Although differing in some details, the exact harmonic method used here is 
essentially the same as a method of Tsypkin (References 75 and 76). 


1. Application of Harmonic Method 

A block diagram of the on-off control system is shown in Figure 6.2, and it is clear that 

£s) = G(s)F(s). (6.3) 


*One-sided limit cycles are possible because the uncontrolled system is unstable. 
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From Equations 6.1 and 6.2, it is found that the 
transfer function is 



Figure 6.1-On-off control characteristic. 


(s + X)\s 2 + (B l - 1)1 

G(s) = ■ 7 =4 4. (6.4) 

For convenience, g(s) = 1/G(s), and Equation 6.3 
becomes 


g(s)£(s) = F(s). (6.5) 

Typical oscillation waveforms for a stable limit 
cycle are depicted in Figure 6.3.* In this figure 

6 - cot, (6.6) 



G(s) 


F(s) 



e<s) 




f i 

k 







Figure 6.2-Block diagram of the on-off control 
system. 


where co is the limit-cycle frequency. The Fourier 
representations of £ and F can be written 


and 


CO 

£(6) = a Q + sin d + ^ ^ sin ( nd + i// n ) (6.7) 
n-2 

b rc 

F(d) =-^-+ ( b n cos nd + c n sin nd). (6.8) 


Using the waveform for F(d) in Figure 6.3, the 
Fourier coefficients of Equation 6.8 are found to be 


if 

J o 

’»=7 / 


K ox( 0 Z - 0 l) 

F(d)dd = — , (6.9) 


b n / F (°) C0S n0 d0 

'0 


K, 


Utt 


(sin n0 2 - sin nd p, 


( 6 . 10 ) 


and 



F(d) sin nd dd - 


cx 

— (cos nd. - cos nd 0 ). 
nn 1 d 


( 6 . 11 ) 


*Notice that there are discontinuities in t whenever the control is turned on or off. Tsypkin’s graphical method as 
given in References 75 and 76 is not applicable to this special case. However, a modified method, which can be ap- 
plied to problems with discontinuities of this type, has been given by Korolev (Reference 77). 
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Writing 


g(s) = g(jco) = R(o))e^^ (6.12) 



and taking the inverse transform of Equation 6.5 
gives 

g(0)a Q + R((z))a 1 sin |# + <£(co)J 

00 

V R(nco) a n sin ^nd + if/ n + = F(d). (6.13) 


jLj 

n=2 


Equation 6.13 leads to the relations 
K c X (°2 ~ 


2ng(0) 


2 K 


cx 


a 


cx 


1 rrR(co) 
2 K 


sm 


e 2 e x 


cx n (@2 

sin — 


o, e 2 

Figure 6.3-Oscillation waveforms (6 = cot). 


n nrrR(nco) 

~ (#i + %)]> 


(6.14) 

(6.15) 

(6.16) 
(6.17) 


and 


=' | + 0 g )J - <j6(acc). (6.18) 

From Equations 6.4 and 6.12, it is found that 


g(jco) 


(2B l + 1 )(B l - 1) - (J5 l - 2)co 2 - 


R(co) = 


and 


(A + /'gj)|(£? l - 1) - &> 2 J 

(6.19) 

- ( 2B l + 1 )(B l - 1) ( w 2 - a 2 )(co 2 + /3 2 ) 

1)]7a 2 + CO 2 (co 2 ~ y 2 )yjx 2 + co 2 

(6.20) 

tan cj)(co) = 

(6.21) 


The constants a and ft can be obtained from Table 2.3. a > y, and it is apparent that R(co) < 0 for the fre- 
quency range y < co < a. This frequency range is not examined in this study, and it will always be assumed 
that R(co) >0.* Equations 6.14, 6.17, and 6.18 can now be written in the form 

In the frequency range y K co < a, the control is turned on when £ is decreasing and is turned off when is in- 
creasing. Therefore, the presence of a stable limit cycle in this frequency range is highly improbable. 


75 



K c^2 ~ W 

3 0 2tt(2B l + 1) ’ 


( 6 . 22 ) 


(d 1 + 0 g ) = 77- + 2 tan 1 (y*)> 

(6.23) 

tn ^ + 0 2>] + tan_1 (x) 

(6.24) 


From Figures 6.1 and 6.3, it is clear that the deadzone and hysteresis constants for the on-off control char- 


acteristic are given by 


and 



(6.25) 

(6.26) 


Finally, with the aid of Equations 6.1 and 6.2, it can be shown that 


m = X + sin [<> - tan_1 (f)] 

Y \ c + (O 


\A 2 + ( n ^f 


sin ^nd + if/ n - tan 1 


(6.27) 


and 


m = 


2a 


cos ^0 - tan 1 (y-) 


p - (B L - i)]v^® + 

■ °° ^ 2fia^co 

4^2 [ c n ") 2 ~ ( b l ~ 


cos ^nd + if/ n - tan 1 (6.28) 


Notice that K cx is just a scale factor, and the “normalized” limit cycle is completely determined by choos- 
ing the parameters 

(0 - 0j), where (0 2 ~ 0^)/2tt = m is the fraction of the total time that the control is on, 

A., where A is the slope of the switching line in the x-x phase plane, and 
co, where co is the limit-cycle frequency. 


2. Closed- Form Solution 

Although the results of the previous section are exact, they are in the form of infinite series. With 
results in this form, the synthesis of the control is sometimes inconvenient. Fortunately, some simplifica- 
tion is possible for the important special case where the control is on for half of the limit-cycle period and 
off during the other half. That is 
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0 2 - 6 1 = 7T, 


(6.29) 


and with Equation 6.16, Equation 6.7 becomes 


2 K„ v (-D” sin [( 2n + D 0 + '/'2m- 1] 

&8) = a 0 + a j sin 8 + — ^ (2n + l)R[(2n + ]>] 

n-1 


(6.30) 


From Equation 6.23 


and 


■(f> 


6 2 - n + tan 


(6.31) 


and 


and 


sin $ i = 


<D 




2 + <o 2 


sin = - 


co 




ar 


Similar expressions can be found for Equation 6.24, and it can be deduced that 


and 


-< 

* 

o 

2 K cx co yy 

[(2/1 + 1 ) 2 co 2 ~ y 2 ] 

8(2 B l + 1) 

17 4^[ (2n 

+ 1 ) 2 co 2 - a 2 ] [(2/i + l) 2 w 2 + /3 2 ] 

_ , K c ^ 

2 K cx a> y^ 

[(2/1 + l) 2 ^ 2 - y 2 ] 

^2 } ~ 2(2 B l + 1) 

17 

+ 1)V - a 2 J [(2n + l) 2 cu 2 + /3 2 ] 


(6.32) 


(6.33) 


Substitution of Equation 6.33 into Equations 6.25 and 6.26 yields 


and 


^ 2(2 B l + 1) 

2K cx a > j a 2 ~y 2 *y 1 

77 (a 2 + ft 2 j(2n + l) 2 - (a/co) 2 ] 


(3 2 + y 2 yy 1 

+ a 2 + /3 2 2Lf [(2n + l) 2 + (j8 /a>) 2 ] 


(6.34) 


(6.35) 


Closed-form expressions for the two infinite series are given in Reference 78, and Equation 6.35 becomes 
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(6.36) 


K 


8 = 


cx 



\a(ar +13*) 


(3(a* + (3*) 


tanh 



A graph of Equation 6.36 for B h = 4 (a 2 - 4.2915, /3 2 = 6.2915, y 2 = 3) is given in Figure 6.4. 

It can be shown that the maximum and minimum values of x(<9) occur at'0 = 0 X + (n/2) and 6=6 + 
(n/2), respectively. Substituting these values into Equation 6.27, it is found that 


and 


where 



y 

K cx „ 

(6.37) 

max 

2(2 B l + l) + a 

X min 

K cx w 

3(2B l + 1 ) «• 

(6.38) 


(-1)" |(2n + l)^ 2 - y 8 ] 

(6.39) 

(2J 1 + 1)[( 

[2n 4- l) 2 co 2 - a 2 ] [(2/2 + 1 ) 2 <o 2 + ft 2 ] 


The summation formulae of Reference 78 can be used to write Equation 6.39 as 


2 K 


cx 


2 _ ,,2 


H) n 


776 .) 


\a~ + 


" (2n + 1)[(2 n + l) 2 - (a/co) 2 ] 


| P 2 + y 8 (-1)" 

a2 + ^ To (2n + 1 >[( 2n + D a + OS/") 8 ] 


M « 2 -y a 
2 |a 2 (a 2 + /3 2 ) 



/3 2 + y 2 

+ /3 2 (a 2 + /3 2 ) L 1 




7y\ 2 « 2 - y 2 

_\ a /y a 2 (a 2 + f3 2 ) 



/3 2 + y 2 
(3 2 (a 2 + (3 2 ) 


sech 



(6.40) 


A graph of Equation 6.40 for B h = 4 is shown in Figure 6.5. 

Equation 6.38 shows that x min < 0 when H > K CX /2(2B L + 1). Intuitive reasoning suggests that 
limit cycles with x min < 0 are unstable. This conjecture will be substantiated by the stability analysis 
given below. 


B. Stability of Limit Cycles 

The results of the previous section cannot be used with confidence until the stability of the predicted 
limit cycle is verified. Tsypkin (Reference 75) has devised a method that reduces this stability investigation 
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Figure 6.4-Hysteresis versus limit-cycle frequency 
(B l =4 , 0 2 - = rr). 



Figure 6.5-E versus limit-cycle frequency (B L - 4, 
0 2 - 6 1 = n). 


to that of an equivalent sampled-data system. A 
brief outline of this method is given below, and 
quantitative results are obtained for the special case 
where d 2 ~ = n. 

1. Tsypkin’s Method 

Consider the perturbed limit-cycle motion 
of the system shown in Figure 6.2. The input to the 
on-off characteristic is 

m = t LC {t) + £(9, (6-4D 

where £ LC (t) corresponds to the assumed limit cycle, 
and £ (t) represents a small deviation. Similarly, 

F(t) = F LC (t ) + F e (t). (6.42) 

The deviation in the output F e (t) is just a sequence 
of pulses, of height K cx and width At n , which occur 
at the beginning or end of the switching instants t n . 
Assuming that <3C deviation can be 

represented by the impulse train 

CO 

F e (t) - K cx ^2 ■ tn) ’ (6 ' 43) 

where S(t) is the Dirac delta function, and following 
Korolev (Reference 77), A t n is expressed in terms of 
the conditions that immediately precede the switch- 
ing instant 

(6.44) 


^ UQ 

At = — 

n \hc(Q\ 


For the special case where the on and off periods are equal (i.e., d 2 6 1 - n), the switching inter- 
vals are T~ T LC / 2 = n/co, and the deviation is given by 


F e (t) 



n=0 


(6.45) 


Consequently, for small deviations, the stability of the limit-cycle motion where 0 2 0 1 - n can be deter- 

mined by checking the stability of the sampled-data system shown in Figure 6.6. 
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Figure 6.6-Sampled-data system for limit-cycle sta- 
bility test ( 6 2 ~ 0 1 = 77 ). 


2. Results for Special Case (0 g - 0 1 - tt) 

The characteristic equation for the sampled- 
data system of Figure 6.6 can be written in terms of 
the z-transform (z ~ e Ts ) as 


1 4 kGpj(z) - 0 , 


where 


and 



G n (s) = - G(s) = 


m 

Q(s y 


(6.46) 


(6.47) 


(6.48) 


Denoting the poles of G N (s) by p n , the z-transform 
of G n (s ) is (Reference 75) 



P(p n ) 

r eP ° T i 

/a 2 - y 2 \ 

z [cos a T + (A/a) sin a r] - 1 

Q'(p„) 

z - e Pj1 ^ 

\a 2 + /3 2 / 

z 2 - 2z cos aT + 1 


+ 



[cosh f3T + (A/,6) sinh /3 t] 
z 2 - 2 z cosh /3T + 1 



(6.49) 


The denominator of Equation 6.47 is found by the following process. First, Equation 6.7 is differentiated 
with respect to time. This gives 


and at 0 = 0, 


4c(0) “ a i (0 cos 0 + ^ ^ & n nco cos (n0 + ^r n ), 
n=2 


and 


cos <9, ~ 


1_ 


4o(»i) 


2K cx a)A 


[(2fi + l) 2 w 2 - y 2 J 


77 Tlo* + ^ 2c ° 2 ~ a2 \ \@ n + V 2( ° 2 + i^ 2 ] 

K cx A ~ y 2 In a \ /3 2 + y 2 (nft 

2 a(a 2 + /3 2 ) \2<w/ /3(a 2 + /3 2 ) 


(6.50) 


(6.51) 


(6.52) 
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and the Fourier series converges to 


(6.53) 

(6.54) 


(6.55) 

and this expression is used in Equation 6.47. Notice that k is independent of K cx . 

Equation 6.46 is a fourth-degree polynomial in z. Because self-oscillations are examined here, 
one root of this characteristic equation is found at z = -1, and the stability of the assumed limit cycle de- 
pends on the locations of the three remaining roots. If these roots are located within the unit circle in the 
z-plane, the limit cycle is stable. Numerical calculations reveal that the limit cycles are stable for A > 0 
and co > co c , where co c is a function of B L only. For three important values of B L , the critical frequencies 
are 


However, there is a discontinuity in i LC {0) at 6= 6 1 (see Figure 6.3), 
the average value 

, 4c(^)+4c(»P 

4,C<^ P = o • 


It is obvious that 


Therefore 


4 C ^P - &c(« l> = «cx- 


4 c («P = 4c^P + 


K 


cx 


2 ’ 


and 


B l = 3.19042 (Earth-Moon L 2 point) co c - 3.102, 

B l = 5.14760 (Earth-Moon point) co c = 3.522, 

B l = 4.00000 o c = 3.287. 


The variations of the roots for B L = 4 are shown in Figures 6.7 and 6.8. 


3. Procedure for General Case (0 2 - tt) 

The characteristic equation for the general case where d 2 - 6 1 i rr is more complicated than Equa- 
tion 6.46. It is shown in Reference 75 that this characteristic equation is 


1 + 


K 


cx 


K 


cx 


1444)1 i4c<q?i 


Gjl( z ) + 


K cx 




— [Gf (z) - G n (z, m)G N (z, -m)] = 0, (6.56) 


where G N (z) is given by Equation 6.49 with T = T Lc = 2 n/co, m = (d 2 ~ 6^/2n, and (References 75 and 79) 


t n 



P(p„) 

Ze m PnT 


N^z, m) N 2 (z,m ) 


z - e PnT 

— z 

D^z) ' D 2 (z) 


4 



n= 1 


p (pp 

e (l~m)p n T 


N^z, 1 -m) 

iV 2 (z, 1-m) 


z - e PnT 


D^z) + 

d 2 (z) 


(6.57) 


(6.58) 
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Figure 6.7-Root variations as a function of A 
(0 2 - e x - 77 , B l - 4, and a> = 4). 


Figure 6.8-Root variations as a function of co 
(0 2 - ~ 77 , B l = 4, and A = 2). 


/« 2 -y 2N 


cos afflT - cos (1 - m)aT sin amT + sin (1 - m)aT]j, 


(6.59) 


N 2 (z, m) = 


a 45 + /S' 


It can also be deduced that 


^ + Vz cosh (3mT - cosh (1 - m){3T + ^ |z sinh (3mT + sinh (1 /n)/3Tj/, (6.60) 


D^z) - z 2 - 2z cos aT + 1, 
DJz ) = z 2 - 2z cosh /3 T + 1. 


4c^P = 4c^i) 


d - £ (fl •> “ (6.64 

Clc'-'V - Slc'-'V 2 

The roots ot the characteristic equation can now be found. Although eight roots are obtained from Equation 
6.56, four of them are extraneous and should be disregarded. These extraneous roots are just solutions of 
the quadratic equations, Dj(z) = 0 and £» s (z) = 0. Another root is found at z = + 1, and the stability of the 
limit cycle is determined by the three remaining roots. 
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C. Solar Sail Control at the Sun-Earth L 1 Point 

To illustrate the use of the preceding formulae, a solar sail control system for a satellite in the vicinity 
of the Sun-Earth L 1 point is considered here. By merely changing the area of the sail in the prescribed 
manner, the on-off control system shown in Figure 6.2 can be obtained. For the Sun-Earth L 1 point, the 
one-sided limit cycle is located between the Sun and the libration point. 

Consider the case where 0 2 - $ l - n. To make use of graphical results, the constant B L1 is approxi- 
mated by B L1 =4. For B L = 4, oj c = 3.287, and a good stability margin can be achieved by choosing a 
limit-cycle frequency of co = 4 (T Lc = 2n/co = 91.31 days). Figure 6.7 indicates that, for o = 4, small devi- 
ations from the nominal limit-cycle motion are adequately damped by choosing A = 2. The choice of x mi 
which specifies K cx , is usually dictated by measurement accuracy limitations. For x min = 100 km, the 
formulae of Section A. 2 of this chapter give 

K cx = 1.613 x I0" 8 g, x max = 343.6 km, 

= 2.965 x 10“ 6 , and 8 = 4.358 x 1CT 6 . 

Equation 3.1 shows that, for a satellite mass of 10 000 kg and = 2, the area of the solar sail is only 
A s = 172.3 m 2 . 

By using a conically shaped solar sail, it may be possible to stabilize the satellite's attitude as well 
as its position. Analyses of the attitude motion for a configuration of this type can be found in References 
80 and 81. 

D. Limit Cycles: Approximate Analysis 

The results of Section A.l of this chapter are not easy to apply when x max and x min are specified and 

the limit-cycle frequency <o is arbitrary. Therefore, an approximate technique that uses x w and x . as 
, max min 

inputs is presented here. 

By neglecting the gyroscopic coupling in Equation 6.1, a closed-form solution for all values of 0 2 - d 1 
can be obtained. Without this coupling, the y-axis motion is simple harmonic. The remaining portion of 
Equation 6.1 is just 


x - (2 B l + l)x = - F(£>. (6.65) 

This equation is only second order and can be conveniently analyzed in the phase plane. However, before 
the phase-plane analysis is presented, the accuracy of this single-axis approximation is checked for the 
special case where d 2 ~ 6 1 = n. 

1. Accuracy of Single- Axis Approximation 

Proceeding in the same manner as in Section A of this chapter, and denoting single-axis quantities 
with primes, it can be shown that 
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g'(s) = 


( 6 . 66 ) 


( 2 B l + 1 ) - s 2 
s + A 

to 2 + (2 B l + 1) a 2 + e 2 

~V5?T^' 

tan (f)'(coi) = -y. 


(6.67) 

( 6 . 68 ) 


Equation 6.68 is identical to Equation 6.21, and many of the equations of Section A of this chapter are valid 
here. With these relations, it is easy to see that 



K c ^ 

2(2B L + l)’ 


(6.69) 


S' =- 


2K 


cx 

mo 


00 

E 

n = 0 


[(2n + l) 2 + (<r/<y) 2 J 



(6.70) 


x‘ 

max 


K 


cx 


2(2 B l + 1) 


(6.71) 


K. 


and 


2 K 


TTCO 


CX 

2 ~ 


L min 2(2B l + 1) 


("l) n 


K 


cx 


2n + l)[(2n + l) 2 + U/co) 2 ~\ 


2e* 


f - seoh (£)] 


(6.72) 


(6.73) 


With the exception of the quantities S' and H', the solutions for the single-axis and coupled systems are 
equivalent. Therefore, a measure of the accuracy of the single-axis approximation can be obtained by com- 
paring these quantities. The percentage errors for the single-axis solutions are shown in Figure 6.9 for 
B L - 4 (c 2 = 9). Although the discrepancy is quite large at the lower frequencies, it is less than 10% in 
both cases when co > 6. 


2. Phase- Plane Method 


The basic equations for the single-axis approximation are 


and 


x - (2 B l + l)x = - F(£) 


(6.74) 


g = x + Ax. 


(6.75) 


F(£) is shown in Figure 6.1. The geometry for a typical limit cycle in the x-x phase plane is illustrated in 
Figure 6.10. From the figure, it is obvious that 
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SWITCHING LINES 



Between points 2 and 3, F(£) = K cx , and a single integration of Equation 6.74 in this region gives 

x 2 = € 2 (x 2 - x 2 )- 2K cx (x - x 3 ) (6.83) 

and 

x 2 = e 2 q 2 x 2 (l - r 2 ) - 2K cx qx 1 ( 1 - r), (6.84) 

where 

x 3 

>1- (6.85) 

*2 

From Equations 6.79 and 6.84, 


e 2 (q 2 r 2 - 1) 
K cx~ 2 q(r- 1) X l* 


( 6 . 86 ) 


With Equation 6.86, the necessary stability condition, K cx > <? 2 x 3 , can be written in terms of q and r as 


(q 2 r 2 - 1) 
2 q 2 r(r - 1) 


(6.87) 


Integration of Equation 6.83 yields 


l 2S ~ f log 


q 2 r(2 - r) - 1 


(q 2 r 2 - 1) - 2 q(r - l)|q + y/q 2 ^ lj 
Combining Equations 6.80 and 6.88 gives the limit-cycle period (T LC = 2n/aj>) 

[q 2 r(2 - r) - lj [q + y/q 2 - 1 


( 6 . 88 ) 


T LC 2 ^12 + ^3) ~ . log 


(q 2 r 2 - 1) - 2 q(r - l)f<? + yjq 2 - 1 


(6.89) 


3. Sample Calculation at the Earth-Moon Lg Point 

As an example of the usefulness of the approximate formulae given above, consider the on-off con- 
trol of a satellite in the vicinity of the Earth-Moon L g point (B L2 = 3.1904, e = 2.7168). By specifying the 
amplitude ratios, Xg/Xj = 2 = q, and x g /x 1 = 2.5 = qr (therefore r = 1.25), Equation 6.89 gives T Lc - 1.407 = 
6.12 days. The fraction of the total time that the control is on can be determined from Equations 6.80 and 
6 . 88 : 


23 n Htrio 

m =-—= 0.4512 . 
1 12 


(6.90) 


This means that the control is turned on continuously for a period of 2.76 days. For thrust periods of this 
duration, an electrical propulsion device would probably be utilized. Measurement accuracy considerations 
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influence the choice of x 1 (x mi ). If Xj = 10 km, then x g = 20 km, x g = 25 km, and Equation 6.86 gives 
K = 2.798 x 10“ 7 g. However, the average control acceleration is only K cx = tnK cx = 1.262 x 10~ 7 g. 
Finally, the slope of the switching line is taken as A = 3, and simple calculations yield = 1.561 x 1CT 4 
and 8= 1.224 x 10" 4 . 

The values listed above for K cx , A, and 8 were used in an analog computer simulation of the 
coupled system (Equation 6.1). It was found that the error in T LC was about 14.1% of the true value 
{T t r) = 1-233 = 5.36 days. Additional results of this simulation are shown in Figure 6.11. Notice 

that x g /x 1 is very close to the value used in the single-axis approximation ( qr = 2.5). 





x versus x. y versus y. 

Grid size is 0.027 m/sec x 10 km. 


y versus x. 

Grid size is 10 km x 10 km 


Figure 6. 11- Analog computer simulation of Equation 6.1 for sample calculation at the Earth- 
Moon L g point. All trajectories are clockwise. 
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CHAPTER VII 


STABLE CABLE 


Consider two satellites of equal mass that are connected by a light cable of adjustable length. Is it 
possible to stabilize the position of the mass center of this configuration in the vicinity of a collinear libra- 
tion point by simply changing the length of the cable with an internal mechanism? This intriguing possibil- 
ity is examined below, and it is found that positional stability can be achieved by varying the length of the 
cable in a manner that takes advantage of the nonlinearities of the modified potential field in the vicinity of 
the libration point. 

This problem was originally analyzed by Colombo (Reference 3), who also concluded that appropriate 
variations in the length of the connecting cable would produce positional stability. However, Colombo’s 
analysis considered only linear variations of the modified potential field in the vicinity of the libration 
point, and it can be shown that his conclusion was based on a fallacious argument. (The adjustable length 
cable was treated as a workless constraint.) 


A. Stabilization Procedure 


In this section, the control logic for the cable variation is formulated, and the stability of the resulting 
system is examined. The geometry for the cable-connected satellite in the vicinity of a collinear libration 
point is shown in Figure 7.1. If the mass of the cable is neglected, the equations of motion for the cable- 
connected satellite are approximately (terms higher than second order are neglected in Equation 2.43, and it 
is assumed that p = v - 0)* 


.x'i - 2 y\ ~ (2 B h + l)x'j ±f-C L [2x'i 2 - (yf + zf)] = F x , 
*2 - 2 ^2 - (2 B l + l)x' 2 ±f C l [2x!, 2 - (y'g + zg)] = - F x , 


and 


7 i + 2x\ + (B l - 1) y\ + 3 C L x' 1 y' 1 - F yl 
y' 2 + 2x(, + (B l ~ 1 )y 2 + 3C l x^ = - F y , 


+ ^L X l Z i ~ F-Z* 


z 2 + B L Z 2 + 3C L X 2 Z 2 “ F z’ 


(7.1) 


*The upper signs hold at Lg, and the lower signs at 
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Figure 7.1-Cable-connected satellite in the vicinity 
of a collinear libration point. 


where 

F = F x i + F y j + F z k, (7.2) 

and the tension in the cable is 

T = mF. (7.3) 

Although fluctuations of the tension will occur when- 
ever the cable is extended or retracted, it is assumed 
here that some tension will always be present; i.e., 
the cable does not become slack. 

1. One-Dimensional Analysis 

To gain some physical feeling for the sta- 
bilization procedure, an approximate, one-dimensional 
analysis is presented before considering the more 
complicated, three-dimensional system. 

Shifting equilibrium point 

Neglecting all coupling terms and as- 
suming that the motion is limited to the x'-axis 
reduces Equation 7.1 to 


*i - (2 B h + l)x'j ± ZC L xf = F x 

and (7.4) 

+ iyx' g ±3C L x'* = -F x . 

Elimination of F x gives 

*1 + x 2 “ ( 2B l + l)( x i + x 2) 1 3 C l( x i S + x 2 2 ) = (7.5) 

Assume for the moment that the mass center of the satellite is at x[, m = ± a, and that the cable length i is a 
constant (f = a). In this case 


and 


x' 



a 



a. 


For equilibrium x^ = x' 2 = 0, and using Equations 7.5 and 7.6, it is found that 


a - 


(2B l + 1) 

" 6C l ' 


1 

2 


(2 B l + l) 2 

9C| 


- a* 


l/S 


3 C, 


4(2 B l + 1) 


(7.6) 


(7.7) 
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Therefore, because of the nonlinear effect, the equilibrium point for the cable-connected satellite is depend- 
ent on the length of the cable. 

To see how stabilization of the mass center might be attained by taking advantage of this con- 
sequence, consider a satellite near an L 2 point with an initial cable length of £(0) = a., and a corresponding 
equilibrium point at x' eq (0) - a t . Because of the natural instability of the instantaneous equilibrium point 
x' the mass center will always experience a repelling acceleration relative to this point. Now suppose 
that the mass center is initially located at x{, m (0) = oq + e(0), and assume that x^ m (0) = HO) = 0. If HO) > 0, 
the initial acceleration is xj, m (0) > 0, and x^ m will increase. However, by extending the cable, x^ q can be 
increased, and it may be possible to shift the equilibrium point to the other side of the mass center; i.e., at 
some time t f , x 1 (t f ) > o i + t(t f ). This shift would decelerate the mass center, and if HO) is not too large, 
the motion of the mass center could be reversed. If e(0) < 0, the initial acceleration is x^ m (0) < 0, and x^ m 
will decrease. With this initial condition, the appropriate shift of the equilibrium point would be achieved 
by retracting the cable [notice that the equilibrium point cannot be shifted to the other side of the mass 
center when e(0) < -oq]. 

Naturally, the magnitudes of the initial errors [e(0), HO)] that can be tolerated with this stabi- 
lization method are bounded. This limitation will be discussed in Section B.l of this chapter. 


Cable control 

For convenience, the coordinate translation 

x' = x ± a (7.8) 

is introduced. With this translation, Equation 7.5 becomes 

x’i + x'g - ( 2B l + l)(x 1 + x g ± 2a) ± 3C l [x^ + x| ± 2a(x 1 + x 2 ) + 2 a 8 ] = 0, (7.9) 


and the equilibrium solution is 


and 



(7.10) 


This equilibrium solution is unstable if the cable length is constant, but the situation is quite different 
when the cable length is varied according to the control law 

i = Xg - Xj_ = [a ± 5(x 1 + Xg) + c(x 1 + Xg)] . (7.11) 

The form of this control law has been influenced by the qualitative arguments of the previous section. 

Stability conditions for the constants a, b, and c can be determined by investigating small 
motions about the equilibrium solution of Equation 7.10. Substitution of 
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and 


cL u 

1 ~ 9 ^ *=1 


a >. 
x p -7T+ 


(7.12) 


2 ^ 

into Equations 7.9 and 7.11 [neglecting terms of 0(<f 2 )] gives 

?! + - (2 B l + l)(f ! + fg) ± 3C L [a(f 8 - fj) ± 8a(fj + fg)] = 0 

and 

(?* - fj) = ± + f 2 ) ± c(?j + ? g ). 

The coordinate for the mass center is 

x l + x 2 £i + ^2 


(7.13) 

(7.14) 


x cm ~ 


2 2 

From Equations 7.13, 7.14, and 7.15, it is clear that 

*cm + 3C L ac * cm +[sC L (ab + 2a) - (2 B h + l)|x cm = 0. 

Therefore, the mass center is asymptotically stable if 

ac > 0 

2B l + 1 


(7.15) 


and 


(a b + 2 a) > 


3 C, 


(7.16) 

(7.17) 

(7.18) 


The constant a can usually be neglected in Equation 7.18, and this stability condition can be approximated 

by 

2B l + 1 


a b > 


3C, 


(7.19) 


For the Earth-Moon L g point (B L2 = 3.1904, C L2 ~ 15.845), Equation 7.19 gives a b > 0.155. 

2. Three-Dimensional Analysis 

Using the coordinate transformation (cf. Equation 7.8) 

x 1 - x ± a, y' - y, z' ~ z 

in Equation 7.1 gives 

x x - 2 y 1 - (2 B l + l)(x 1 ± a) ±-|c l |2x| ± 4 ax 1 + 2a 2 - ^y 2 + zfj^ = F x , 


(7.20) 
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- 2y, - (ZB, + l)(x 


0 ±l c f-N 


4 aXg + 2a ^ 


-(y! + zi)] = 


y\ + 2x : + (B l - l)y 1 + 3 C L (x 1 t a)y 1 = F yt 

y 2 + 2 x g . + (B l - l)y 2 + 3C L (x 2 ± a ^2 = “ F y’ 

z 1 + B l z 1 + 3C L (x t ± a)z 1 = F z , 
x 2 + S L Z 2 T 3C l/ X 2 1 a ^ Z 2 = " F z’ 


where a is given by Equation 7.7. For a cable length t = a, it is readily verified that an equilibrium sola 
tion of Equation 7.21 is 


and 


a a_ 

x l~"2’ x 2~2 

y l = y z = z l = Z 2 =0. 


(7.22) 


This equilibrium situation is depicted in Figure 7.2. Notice that stabilization about this equilibrium solu- 
tion must control the relative orientation of the cable-connected satellite as well as the position of its 


center. 

The cable control law for the three-dimensional case is not obvious, but after a few trials, a satis 
factory control law was devised. This control law is given by 


e 8 = (x 2 " x i) 2 + " ^l) 2 + ( z 2 ~ z l f 

- ja ± b(x 1 + x 2 ) ± c(x 1 + x 2 ) + /3(y ^ - y g )J > 


where a, b, c, and j8 are constants. It is hoped that the reasoning behind this particular formulation will be 
made clear by the stability investigation presented below. 

Before proceeding to the stability analysis, some useful equations are obtained. Using Figure 7.1, 
it can be verified that (check Equation 7.2 and 7.20) 


F x X 2- X l 

(7.24) 

F y-y z -y\' 


F x X 2' X 1 

(7.25) 

F z~ Z 2~ Z l' 


I ^ 

1 

ro 

I 

'c 

t-*- 

(7.26) 

F z z 2 ~ Z 1 



Rearrangements of Equation 7.21 give 
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(7.27) 


X 1 + x 2 2(y 1 + y g ) - (2S L + l)(Xj + Xg ± 2a) ±-f C L [2^ + x|) ± 4 a(x 1 + x g ) 

+ 4 a 2 - (y^ + y 2 + z^ + z|^]= 0 , 

Jl + y g + 2(x a + x g ) + (B l - X)(y 1 + y 2 ) + 3 C L [(x 1 ± a)y 1 + (x g ± a)y g J = 0. (7.28) 

Z 1 + z 2 + B l ( z i + z 2 ) + 3 C l[( x i ± «) z i + ( x 2 ± a)z g l = 0, (7.29) 

- x‘i - x‘g - 2(yi - y g ) - (2 B l + l)(x t - x g ) ±-gC L j 2 ^x 2 - x|^ 

± 4a( Xj - Xg) - [(y 2 - y 2 ) + (z\ - z|)]|, (7.30) 

2F y = 4 ~ 4 + 2 <4 - * 2 > + ( B l - lXy 4 - y s ) + 3 C l [(Xj ± a)y t - (Xg ± a)y g J, ( 7 . 31 ) 

and 

2F Z = z t - z 2 + B l (z 1 - z 2 ) + 3 C l J(x 1 ± «)z 1 - (x g ± a)z 2 J. (7.32) 

For small motions about the equilibrium solution of Equation 7.22 


a , ^ a ^ 

x l _ 2 + ^l’ X 2~Y + ^ 2 ’ 


yi = 7r ^ (7.33) 

Z 1 = ^1’ Z 2 = '=2' 


Substituting Equation 7.33 into Equations 7.23, 7.27, 7. 
y, 0 gives 



Figure 7.2-Equilibrium situation for the cable- 
connected satellite near an L g point. 
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5, and 7.29, and retaining only linear terms in (£ 

i - 4> ± *4 + * «4 + 4) 

+ /3(,j - , 2 ) = 0, (7.34) 

. + 4 - 2«i + ijg) - (2 B l + l)(f x + fg ± 2a) 

± l C t[ a8+2a( f 2 “4) 

± 4a(fj + <f 2 ) + 4a 2 J = 0, (7.35) 
+ i/g + 2(4 + 4) + (B t - l)(i), + r/g) 

+ 3C l ^-|-± a^j +^| ± ajijg = 0, (7.36) 



and 


?! + ? 8 + b lUi + Q + 3 C l[(“ f ±a k + (l 1 a )^] = °- 


(7.37) 


Two more equations are obtained by eliminating F x , F y , and F z in Equations 7.24, 7.25, 7.30, 7.31, and 
7.32. Again keeping only linear terms in (£ rj, 0> this yields 


and 


rj 2 ~ r} 1 + 2(4 - +|' c L a(iji + ^ + 3C L a )^2 

? g - ^ Tfc^ + Cg) + [8(B t - 3C L a) + l](f B - = 0. 


(7.38) 

(7.39) 


Equations 7.34 to 7.39 constitute six independent equations for the six unknowns (£ r £ 2 , r/ 2 , ^ £ 2 )- 

These equations will be used to test the stability of the controlled motion of the cable-connected satellite 

about the equilibrium solution of Equation 7.22. 

The linearized motion out of the xy-plane is uncoupled, and the characteristic equation for this 
motion can be determined from Equations 7.37 and 7.39. A simple calculation gives 

s 4 + a^s 2 + a g = 0, (7.40) 


where 

= (4 B h + 1) - 12C L a (7 - 41) 

a nd , 

a 2 "f C l/ 2a 4 a )][( 3B L + "f C l/ 6a 4 a )J 

+ [B L ~§C l (2« ± a)][(3B t + 1) -|c t (6a + a)]| . (7.42) 

Terms involving a and a 2 can usually be neglected, and the roots of Equation 7.40 are approximately ± iV&L, 
± j\/3 B l + 1 (j = /I). Therefore, according to the linearized analysis, the out-of- plane motion of the cable- 
connected satellite is bounded. However, small oscillations, which could be caused by initial-condition 
errors, cannot be damped by the cable control. 


The characteristic equation for the coupled motion in the xy-plane is obtained from Equations 
7.34, 7.35, 7.36, and 7.38. After some tiresome manipulations, it is found that 


where 


6 + b t s b + bgS 4 + h 3 s 3 + 5 4 s 2 + b 5 s + b 6 = 0, 

(7.43) 

b x ~ 3C L ac - 2/3 , 

(7.44) 

b 2 = (k 2 + k 3 + 4 - kj + 3C L ab, 

(7.45) 

b 3 = 3 C L (k 2 + k s + 2)a c + 2^ - (k g + 4)]/3, 

(7.46) 
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(7.47) 


and 


and 


and 


b 4 = k &( k 2 + 4 ) ~ k i ( k 2 + k 3 > ~\ C L a2 +'3C L (k 2 + k s + 2 )ab, 


^5 3C L (k 2 k s 4 C^a^jac + 2^k 1 k 2 + g 

(7.48) 

b 6 = (^2 k 3 C£a^(3C L ab - kj, 

(7.49) 

k \ - (2 #l + 1) “ 6 C L a, 

(7.50) 

k 2 ~ ( B L “ 1) “ 3C L a, 

(7.51) 

kg = 3B^ - 9 Cj^a. 

(7.52) 


The determination of suitable values for the control parameters b, c, and /3 can be accomplished in an effi- 
cient manner by employing the root-locus technique. Taking c = /3 = 0, Equation 7.43 is written in the form 


where 


and 


3C L ab|s 4 + c^s 2 + cj 


s 6 + djS 4 + d g s 2 + dg 

(7.53) 

= k 2 + kg + 2, 

(7.54) 

- k k 

~ K 2 K 3 4 L 'L a ’ 

(7.55) 

= k 2 + k s + 4 ~ k v 

(7.56) 

= - jk x (k 2 + k g ) - kg(k g + 4) - -C% a 2 ], 

(7.57) 

“ “ k l(^2 k 3 ~ 4 C L a2 )' 

(7.58) 


A root-locus plot of Equation 7.53 for the Earth-Moon Lg point is given in Figure 7.3.* The root-locus plot 
shows that the system is neutrally stable when ab > (ab) 3 - 0.600. Notice that (ab) 2 corresponds to the 
critical gain for the one-dimensional case, i.e., (a b) 2 = (2B L + l)/3 C L (cf. Equation 7.19). The critical 
gains (ab) g and (a b) 3 are listed in Table 7.1 for several libration points. In every case, (ab) 3 /(ab) 2 < 4, and 
it can be concluded that neutral stability is insured at any libration point when 


It was found that the contributions of terms involving a and <x were negligible for the values of a that are used in 
the present study With a = 0 in Equations 7.54 to 7.58, the critical gains are (a b) 1 = 0.084, (a b) 2 = 0.155, and (ab) 3 = 
0.600. For a = 10 in Equations 7.54 to 7.58, the critical gains are (ab) 1 = 0.085, (ab) 2 = 0.155, and (ab) 3 ~ 0.599. 
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a b > 


(7.59) 



Figure 7.3-Root-locus plot for Equation 7.53 at the 
Earth-Moon Lg point (B L2 = 3.1904, C Lg = 15.845). 
Gain is a b. 


and 


4 (2B l + 1) 
3Cr 


Damping is introduced by choosing appro- 
priate values for the control parameters c and (3. 
For this case, Equation 7.43 is written as 


where 


and 


3C L ac^b' 1 s 4 + 5gS 2 + b'^s 
s 6 + 5gS 4 + b 4 s 2 + 5 6 

b\ = 1 - 20', 

hg = (/fg + /fg + 2) 

+ 2 \k t - (. kg + 4)]j8', 

+ 2 ( k i k z +l c E a2 V 


"I, 


(7.60) 


(7.61) 


(7.62) 


(7.63) 


0' 


i9 

3C L ac ‘ 


(7.64) 


It should be emphasized that both c and /3 are needed in Equation 7.43 to obtain a damped system. The 
parameter /3' (which determines /3) must be chosen so that the zeros in Equation 7.60 are located between 
the poles on the imaginary axis. (This pole-zero configuration is not possible when /3' = 0.) If this is not 
done, a segment of the root locus of Equation 7.60 will be in the right half of the s-plane. 


Table 7.1— Critical gains for cable-connected satellite. 


Libration point 

(ab) 3 

(ab)g 

(ab) 3 

(ab) 2 

Earth-Moon L 1 

0.651 

0.175 

3.72 

Earth-Moon L 2 

.600 

.155 

3.87 

Sun-Mercury L g 

.0144 

.00378 

3.81 

Sun-Venus L 2 

.0355 

.00932 

3.81 

Sun-Earth L g 

.0380 

.00999 

3.80 

Sun-Mars L g 

.0181 

.00475 

3.81 

Jupiter-Io L g 

.0882 

.0232 

3.80 

Jupiter-Ganymede Lg 

.114 

.0298 

3.83 

Jupiter-Enceladus L 2 

.0140 

.00369 

3.79 
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By fixing the values of ab and (3 1 , a root-locus plot of Equation 7.60 for the gain a c can be con- 
structed. The plot for the Earth-Moon h 2 point with a b = 0.750 and /3' = 0.200 is shown in Figure 7.4. 
Although more damping would be desirable, the system is asymptotically stable for all values of the gain a c. 
When ac - 0.100, the roots are located at -0.0475 ± j 1.2249, -0.7595 ± j 3.8841, and -0.4764 ± 74.9958, and 
from Equation 7.64, /3 = 0.951. 

B. Other Considerations 

Although the stability of small motions about the equilibrium solution of Equation 7.22 has been estab- 
lished, several questions concerning the theoretical feasibility of the cable control system are still unan- 
swered. Some of these questions are considered in this section. 


1. Cable Extension Limitation 

The extension of the cable is made possible by the gravitational gradient in the vicinity of the 
libration point. Therefore, the acceleration of the cable’s length i is bounded [£ max = (2 B h + 1)1]* This 
restriction has the effect of placing upper bounds on the initial errors that can be tolerated with the cable 
stabilization technique. In this section, approximate analytical estimates of these error bounds are obtained 
by using the one-dimensional model of Section A.l of this chapter. 


tot 


Denoting initial values with a subscript i, the coordinate translation of Equation 7.8 is specialized 

x' = x + aj . (7-65) 


The coordinates of the two masses are 


and 


Xj_ = - q + e 


(7.66) 


x 2 = q + 6. 

Initially, it is assumed that q(0) = a i /2, <?(0) = 0, and e(0) ~ o( af). Substitution of Equation 7.66 into Equa- 
tion 7.4 gives (using Equation 7.65) 

- q + V - (2 B L2 + l)( f - q + a t ) + 3C L2 [q 2 - 2 qe + f 2 + 2a ^ - q) + af] = F x (7.67) 

and 

q + e - (2 B L2 + 1 )(q + e + a t ) + 3C L2 jq 2 + 2 qe + <r 2 + 2 a^q + e) + af] = ~F X - (7.68) 

If terms of O^af^ are neglected, Equation 7.67 or Equation 7.68 yields 


q-(2B L2 + 1 )q = -F x . 


(7.69) 


*When the cable is retracted, ^ max is determined by the tensile strength of the cable. 
TFor clarity, double signs are deleted, and the analysis is limited to an Lg point. 
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Figure 7.4-Root-locus plot for Equation 7.60 at the 
Earth-Moon L 2 point (B L2 = 3.1904, C L2 = 15.845). 
Gain is a c, a b = 0.750, and j8' = 0.200. 


During the cable extension period (0 < t < t^) it is 
assumed that 

F x = (X- y)(2B L2 + 1)9 = 0, (7.70) 

where 0 < y < 1. Equation 7.69 can now be written 
as 

q - y(2B L2 + 1)9 = 0, (7.71) 

and with the initial conditions given above, the 
solution is 

a i 

q(t) =-g-cosh kt (7.72) 

and 

ka i 

q(t)=—sinhkt, (7.73) 

where k 2 = y(2B L2 + 1), When the cable is fully ex- 
tended, q(tj) = a m /2, and the extension time t 1 can 
be obtained from 


cosh kt 1 =jr, (7.74) 

where 8 = a, i /a m (0 < <5 < 1). 

To find the motion of the mass center, add Equations 7.67 and 7.68 and neglect terms of o(af). 
This yields (using Equation 7.7) 


€ - (2B ^ 2 + l)f + 3C^g 




(7.75) 


For the cable extension period (0 < t < tj), Equation 7.75 is 

- (2B L2 + l)e + ~|c L2 § 2 a 2 (cosh 2 kt - 1) = 0. 

The solution of Equation 7.76 is given by 

e(t) = [^(O) + (A + J5)] cosh k't + -~r sinh k't - (A + B cosh 2 kt), 

where 

k' 2 ^(2B L2 + 1 ), 


(7.76) 


(7.77) 
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A = 


3C L2 S8a m 

8(2£? L g + 1) ’ 


and 


B =■ 


3 C L2 S^ 2 m 


8(SB L2 + l)(4y-l)- 
After the cable has been fully extended (t > tj), Equation 7.75 becomes 

v “ ( 2B L2 + +_ | C L2 a ffl (1 " = °’ 


The solution of Equation 7.78 is 

T 1 4^) . 

e(£) = Utp - A'\ cosh k'(t - tj + —p~ sinh k (t - t x ) + A , 


where 


A' = 




4(2B L2 + 1) 

The maximum allowable initial conditions [e(0), e(0)] are determined by assuming that, 


time t c 


and 


e(t 2 ) = 0 


<t 2 ) = a m - a t = 


3C L2 3 m^ 1 ^ _ A , 

4(2B L 2 + 1) 


= A 


Using Equation 7.80 in Equation 7.79, it is readily deduced that 

4^) 

tanhk'(tg-t 1 )= F prr^ 


= 1. 


Finally, with the aid of Equations 7.74, 7.77, and 7.81, it is found that 


40) + 


40) 


3C L2 3 m 


V2B Lg + l _ 8(2B L 2 + l) 


G(y, 5), 


where 


G(y, 8) = e k |(2 - 8 2 ) + j^(cosh 2 kt ± + 2 V7 sinh 2 fc£ i)] 


4 y§' : 


(4y - 1) 


When f(0) = 0, the ratio <<0 )/c(£ g ) may also be of some interest. From Equations 7.80 and 7.82, 
written as 


(7.78) 


(7.79) 


at some 


(7.80) 


(7.81) 


(7.82) 


(7.83) 
this can be 
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(7.84) 


< 0 ) _ G(y, 8 ) 

«(t 2 ) 2(1 - S 2 ) ‘ 

Graphs of G(y, 8) and «(0)/«(t g ) are given in Figures 7.5 and 7.6, respectively. 

From Equation 7.82 and Figure 7.5, it can be seen that the allowable initial conditions are maxi- 
mized by choosing y = 1, 8 = 0.55, and letting a m become very large. However, other factors call for small 
values of a m and y and larger values of 8 .* These tradeoffs are not investigated here. 

As an example, consider a cable-connected satellite at the Earth-Moon L g point. For y = 0.50, 

8 - 0.8, and a m = 10“ 2 = 3844.05 km, it is found that G(y, 8) = 0.375, a i = 3075.24 km, and a i - 39.61 km. If 
<r(0) = 0, <<0) = 11.60 km, and <t g ) = 22.28 km. If e(0) = 0, e(0) = 0.0839 m/sec. Although this cable is rather 
long, significant reductions of the cable’s length will require high measurement accuracy; e.g., e(0) = 0.725 
km for a m = 961.01 km. 


2. Structural Comments 


The structural aspects of extremely long, cable-connected satellites have received some attention 


in the literature (References 82 and 83). In Reference 
the length of the cable has also been considered. 
These studies have indicated that similar design 
problems for the cable-connected satellite proposed 
in this chapter would not be insurmountable. 

To provide some idea of the magnitudes 
of the structural parameters for a cable-connected 



83, the mechanical design of a device that would vary 



a 


Figure 7.5-The function G(y, 8). 


Figure 7.6-Mass center displacement ratio £(0)/c(t g ) 

U(0) = k(t 2 ) = oi. 


*One of these factors is the impulsive deceleration that takes place when the cable reached its maximum length at 
t = tj_. From Equations 7.73 and 7.74, this impulsive deceleration is given by Axg^i) = 2 a mb / (2B L 2 + “ 8 ^ * In 

actual practice, the deceleration will be distributed over a finite time interval that is determined by the maximum tensile 
strength of the cable. This deceleration period is neglected in the analysis presented in this section. 
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satellite at the Earth-Moon L g point, a few representative calculations are presented here. For the purposes 
of this calculation, it is assumed that 


(1) m = 10000 kg. 

(2) The cable is an aluminum wire (modulus of elasticity = E A1 - 1.03 x io 7 psi, density = p A1 = 
2.70 gm/cm 3 , and ultimate tensile strength = a A1 = 30000 psi). 

(3) The nominal length of the cable is t - 3844.05 km. 

(4) The cable tension is never greater than 10 times its nominal value. 

With these assumptions, a few elementary calculations give* 


and 


£ _ 

T nom = m (^^L2 + ^2~ x d y nes (nominal cable tension), 
O’Al 

a des = = 3000 psi (design stress), 

d Q = 6.22 x 10“ 3 cm (cable diameter), 
m c = 504.15 kg (mass of cable). 


*It has been suggested to the author that the time delay, due to the finite propagation speed of a longitudinal wave 
in the cable, could be troublesome. However, the propagation speed is just y/E A1 /p A1 - 5.12 km/sec, which corresponds 
to a time delay of t ^ = 12.5 minutes for a cable length of £ = 3844.05 km. This time delay is clearly negligible. 
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CHAPTER VIII 


APPLICATIONS 


In this chapter, a number of possible applications for libration-point satellites are presented. Some of 
these proposals are original, while others are derived from previous suggestions. It is not the author’s in- 
tention to present an all-inclusive list, but rather to give some indication of the usefulness of libration-point 
satellites. 

A. Supporting Role for Lunar and Planetary Missions 

1. Utilization of the Earth-Moon Collinear Points in Future Lunar Operations 
Lunar communications 

In the post- Apollo period, with the advent of semipermanent lunar bases and far-ranging sur- 
face vehicles, a capability for real-time communications between widely separated lunar terminals will be 
needed. For bases located on the far side of the Moon, an uninterrupted communications link with the Earth 
would also be desirable. In this section, it will be shown how these communications requirements can be 
satisfied by stationing relay satellites in the vicinity of the collinear libration points of the Earth-Moon 
system. 


The possible use of libration-point satellites for communication between points on the lunar 
surface was first mentioned by Arthur C. Clarke (Reference 9) as early as 1950. However, it has been only 
recently that methods for using a single libration-point satellite to establish a communications link between 
the Earth and the far side of the Moon have been devised. Two of these methods, the “Lissajou orbit” and 
“halo orbit” concepts, were originally proposed by the present author in Reference 7. Another technique, 
the “hummingbird” concept, has been presented by Vonbun (Reference 18). A comparison of these methods 
is given below. 


Lissajou orbit concept. Consider a satellite that is following a quasi-periodic orbit about 
the Earth-Moon L g point. Neglecting higher order corrections, the equations of the orbit are (see Chapter III) 


and 


x n = A xi sin V 


y n = A yl cos V- 


( 8 . 1 ) 
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where A xl = kA yV k = 0.343336, and co n = 1.86265. The geometry for this orbit is given in Figure 8.1. If 
the satellite is simultaneously executing an out-of-plane oscillation, then (cf. Equation 3.29) 

z n = A zl cos (co z t + a Q ), (8.2) 

where co z = 1.78618 and « 0 is the initial phase angle. The trajectory of the satellite as seen from the Earth 
is just a Lissajou curve, and is depicted in Figure 8.2. Because the frequency difference between the peri- 
odic orbit and the z-axis oscillation is very small, the trajectory in Figure 8.2 can be viewed as a slowly 
changing elliptical path. 

To set up an efficient communications link between the Earth and the far side of the Moon, 
the satellite should perform oscillations that are large enough to be visible from any point on the Earth fac- 
ing the Moon (see Figure 8.3). Neglecting the small correction for A xl , the radius of the occulted zone for 
Earth visibility is 


a m^~ y 



Figure 8.1-Geometry in the Moon’s orbital plane. 
Not to scale. 



Figure 8.2-Satellite trajectory about the L g point as 
seen from the Earth. Not to scale. 


+ 3099 km, (8.3) 

where R® = 6371 km and 1738 km (Reference 
28). Slightly larger oscillations by the satellite 
near the L g point will enable it to obtain line-of- 
sight contact with the L 1 point as well as the Earth. 
By augmenting the oscillating satellite with a sec- 
ond satellite stationed at the L ^ point, a surface-to- 
surface communications link between the near and 
far sides of the Moon could be established. (This 
possibility was originally discussed by the present 
author in Reference 16.) Using Figure 8.1, it is 
easy to see that the radius of the occulted zone for 
this case is 

y l, i y y t p 

A ML1 = — Z R D = 3671 km - ( 8 - 4 ) 

y li 

Unfortunately, the oscillating sat- 
ellite will still enter the occulted zone at periodic 
intervals. The fraction of the total time that the 
satellite will be hidden is approximately (using 
Equations 8.1 and 8.2 and Figure 8.2) 



For A m = 3099 km and A yl = A zl = 0.02 = 7688 km, 
this fraction is only r yz = 0.070, i.e., 7% of total 
time. However, occultation periods as long as 
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t y = (2/co a ) sin 1 (A m /A yl ) = 1.94 days will occur 
periodically. Because this time is rather long, some 
means for completely eliminating this occultation 
would be useful. 

The satellite’s motion relative to 
the quasi-periodic orbit can be controlled by any of 
the methods discussed in Chapters IV to VI. If the 
radial-axis control of Figure 5.1 is employed, the 
station-keeping cost could be as low as | F cx | = 

8.86 x 10" 9 g (see Section A.l of Chapter V). Of 
order corrections are included in the nominal path 


Phase-jump control is one method to prevent occultation of the satellite oscillating about 
r Whenever the satellite is about to enter the occulted zone, a control pulse is applied to alter the phase 
angle, but not the amplitude, of the z-axis oscillation so that the portion of the trajectory that passes through 
the occulted zone is bypassed. This change in phase angle is accomplished by impulsively reversing the 
direction of z without changing its magnitude (see Figure 8.4). From Equations 8.1 and 8.2, the parametric 
equations for the trajectory shown in Figure 8.2 can be written as 

y n = A yl cos co n t 
and 

z n ~ A zi cos ("n* + a 

where a = an - et and e = o>„ ~ co_ — 0.07647. Since f is small, the phase angle a will be approximated by an 

U n z 

average value through one cycle. If a c is the phase angle for the cycle that just misses the occulted zone, 
the magnitude of the control impulse is (see Figure 8.4) 

|Az n | =2co n A zl sin a c , (8.7) 

and since pulses must be imparted at intervals of 



Figure 8.3-Geometry for Earth visibility. Not to 
scale. 

course, this low cost can only be realized when higher 
calculation (see Chapter III). 


7 7 “ 2 «„ 

At = 

e 


( 8 . 8 ) 


the average control acceleration is given by 


CZ ' 


2((o n A z ^ sin a c 
n - 2 a„ 


It follows from Equation 8.6 and the geometry of Figure 8.2 that 

2r 2 = (Ay X + A^ + (a% x + A^ t cos 2a) cos 2c o n t - ( A ^ sin 2a) sin 2 (o n t 


(8.9) 


( 8 . 10 ) 
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where r is the radial distance from the origin to the satellite. For a trajectory that just touches the occulted 
zone 


pA £ 


= ( A yl + A l i)~[( 


{A y l + 


A zl cos 2 «c> 


) 2 + sin 2a c ) s ] 


1/2 



1/2 

A 2 - 1 

_A 2 sin 2 a c ~ 1_ 


( 8 . 11 ) 

( 8 . 12 ) 


where A = (A yi /A M ). Substitution of Equation 8.12 into Equation 8.9 yields 


\ F cz\ ~ % €C °n A M^> a c )> 

where 


(8.13) 


K A- a c ) = 


(A 2 - l)^ sin « c 
(77 - 2a c XA 2 sin 2 a c - l)^ 2 


(8.14) 


For a fixed value of A, there exists a unique value of a c , which minimizes the average control acceleration. 
An elementary analysis shows that this optimum value must satisfy the transcendental equation 

2 (a 2 sin 2 a c - l) tan a c = n - 2a c . (8.15) 


The phase angle satisfying Equation 8.15 is given in Figure 8.5 as a function of A. 

Up to now, the higher order corrections to the trajectory of Equation 8.6 have been ne- 
glected. When these corrections are superimposed on the approximate trajectory, the occultation geometry 
will become more complicated. A simple way (possibly not the most efficient) to account for this complica- 
tion is to define a new value of A M (denoted by A‘ M ) in such a way that the true trajectory will not enter the 



original occulted zone so long as the approximate tra- 
jectory does not enter an occulted zone of radius A' M . 



Figure 8. 4- Phase plane representation of the control 
impulse. 


Figure 8.5-Optimum phase angle. 
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In Chapter III, it was found that the amplitude corrections caused by the lunar eccentric- 
ity were much larger than those from other effects. Therefore, amplitude corrections from effects other 
than the lunar eccentricity will be neglected here. Using the results of Chapter III, A M is taken as* 


A M ~ A M + P A y 1 l - 


M 


p\" 


(8.16) 


where A' ^ (A yl /A' M ) and /3 £ 0.17. Replacing A u in Equation 8.13 by A' M from Equation 8.16 gives 

f(A',a c ) 

\F CZ \ - 2eco n A M j _ p-, = 2 € 0 ) n A M g(\ ,a c , (3). 


(8.17) 


The function g(A\ a c , /3) for optimum values of a c is plotted in Figure 8.6 as a function of A . For /3 0.17, 

a nearly optimum value of g( A', a c , (3) is obtained by taking A' = 2. With A m = 3099 km, it is found that 


A‘ m = 4695 km, 


= 43.17°, 


and 


A yl = 9391 km, A zl = 8708 km, 

At = 93.04 days, |Az fl | = 59.08 m/sec, 
l*Q = V-50 x 10"' 7 g- 


Occultation of the satellite oscillating about L g can also be averted 
In this scheme, a single-axis control is used to synchronize the 
In other words, if higher order trajectory corrections are ne- 
glected, the control will produce a closed ellipti- 


Halo orbit concept. 

by using a frequency control technique, 
fundamental y-axis and z-axis oscillations. 



Figure 8.6-The function g( A\ a c , /3) for optimum 
values of a Q . 


cal path in the yz-plane that always avoids the 
occulted zone. Although frequency control can 
be applied in any axis, the cost will be different 
for each case. Minimum costs for three important 
cases are given below. As in the previous sec- 
tion, the effect of the eccentricity correction is 
included by using A' M for the radius of the oc- 
culted zone. 

The frequency control method is 
most easily understood by examining the case where 
z-axis control is employed. In this case, it is as- 
sumed that the approximate satellite orbit in the 
xy-plane is given by Equation 8.1 with A yl > A' M . 
The z-axis control is used to force the satellite to 
follow the path 


*The eccentricity correction for the z-axis oscillation is given in Appendix B. For this case, A M A ^ + /3A zi , 
where /3 = 0. 16. 
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Z n ~ A m COS 


(8.18) 


(i.e., a feedback control is used to zero the coordinate C where z - z n ). The cost for the frequency 
change is given by |F CZ |, where F cz - K cz cos a> n t* Because 

z n + B L2 Z n ~ B cz ~ ^ cz cos (8.19) 


it is easy to show that 

K cz ~ (^n ~ B L2j A M = ("n ~ B L2)~{ - ft 
(/3 = 0.16 for the z-axis). Therefore, the average control acceleration is 




= 0.2115 A„. 


( 8 . 20 ) 


( 8 . 21 ) 


3.10 to 3.12) 


Similarly, for a y-axis control, F cy = K cy cos co z t, and it is found that (using Equations 


[F cy [ = 0.1589 A m . 


( 8 . 22 ) 


With y-axis frequency control, \A xl /A yl \ = 0.3379 (\A X ^/A y ^\ = 0.3433 for the natural oscillatory mode). 

For an x-axis control, F cx = K cx cos co z t, and a simple calculation gives 

|F^| = 0.4703 A m , (8.23) 


with \A xl /A yl \ = 0.2799. 

Taking A Me - 3099 km, Equations 8.21 to 8.23 yield 

|F^| = 4.73 x 10 ~ 7 g, 

\F~^y\ - 3.56 x 10~ 7 g, 

and 

| F cx | = 1.05 x 10" 6 g. 

Hummingbird concept. In both of the concepts presented above, an “orbiting” satellite 
was used to set up a communications link between the Earth and the far side of the Moon. However, this 
communications link can also be obtained by using a satellite that is permanently displaced to one side of 
the Earth-Moon L 2 point. In this case, a continuous thrust control is needed to station the communications 
satellite at the nonequilibrium position. The control acceleration can be minimized by displacing the satel- 
lite along the y-axis, and is approximately (see Equation 5.22 and Figure 5.4) 


*It is assumed here that noise inputs are negligible (see Section A. 2 of Chapter V). 
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(8.24) 


|F c i-|F cy |^(B L2 -l)S y - 2.1904 8 y . 

When 8 y = A m = 3099 km, Equation 8.24 gives |F C | = 4.90 x 10 6 g. 

Although the displaced satellite would be stationary with respect to the libration point, 
its distance to the Moon would vary periodically because of lunar eccentricity. Of course, a constant dis- 
tance between the satellite and the Moon could be maintained by simply canceling the eccentricity effect 
with additional thrust control. A rough calculation shows that the average control acceleration for this task 
would be about |F C | = 1.20 x 10 1 2 3 4 5 g . 

Comparison of different techniques. It is clear that, of the three procedures considered 
above, the halo orbit concept with y~axis control is the most economical method for maintaining a continuous 
Earth-to-lunar far-side communications link. However, the cost for the halo orbit concept with z-axis con- 
trol is only slightly higher, and other factors may influence the choice between the y-axis and z-axis con- 
trols. The cost for the halo orbit concept, with either y-axis or z-axis control, is lower than the cost for the 
hummingbird concept by an order of magnitude. 

If occasional occultation periods can be tolerated, the Lissajou orbit concept should be 
considered. Uninterrupted communications would still be possible for intervals of about 3 months, and the 
station-keeping cost could be reduced by an order of magnitude. (This is two orders of magnitude below the 
cost for the hummingbird concept.) For emergency situations, when longer continuous communication inter- 
vals are required, the phase- jump control technique could be employed. 

Concluding remarks. Methods for using a single libration-point satellite to obtain con- 
tinuous communications between the Earth and the far side of the Moon have been presented. The possibil- 
ity of using two libration-point satellites (one oscillating about the Lg point while the other is stationed at 
L 1 ) to establish a point-to-point communication link covering most of the lunar surface has also been men- 
tioned. The libration-point communication system would be useful during the early lunar exploration period, 
but its full potential will probably not be realized until more advanced lunar surface operations are carried 
out. For these more advanced missions, the libration-point network could also function as a navigation and 
control center for manned and unmanned surface vehicles. 

Other proposals (References 84 to 86) for obtaining Earth-to-lunar far-side and long-range 
lunar surface-to-surface communication have recommended the use of relay satellites in lunar orbit. How- 
ever, orbiting relay satellites have several obvious disadvantages. Some of these disadvantages are 

(1) Many satellites are needed for adequate coverage. 

(2) Tracking and acquisition problems are difficult. (These problems are considerably 
easier with libration-point satellites, because these satellites are almost stationary with respect to the 
lunar surface.) 

(3) The contact time for any given satellite is relatively short, and frequent switchovers • 

are necessary. 

(4) The diverse antenna pointing requirements (satellite to Earth, satellite to lunar sur- 

face, and satellite to satellite) would complicate the attitude control problem. 
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Although orbiting satellites possess some advantages over libration-point satellites (e.g., communication 
distances are shorter and station keeping is usually not required), they do not appear to outweigh the disad- 
vantages listed above. A recent investigation by Schmid (Reference 87) has also expressed a preference for 
the libration-point satellite technique. In Schmid’s study, a quantitative comparison of lunar far-side com- 
munication requirements (power, beamwidth, antenna gain, bit rate, etc.) for a libration-point satellite and a 
1000-km altitude satellite in a circular lunar orbit is presented. 

Rendezvous technique 

Present plans for the Apollo lunar landing mission call for the insertion of the combined 
spacecraft [command module (CM) and lunar landing vehicle (LLV)] into a low lunar orbit. The LLV will 
then descend to the lunar surface, stay there for 1 to 2 days, ascend to the lunar parking orbit, and rendez- 
vous with the CM. The CM will then return to the Earth. 

An alternate method uses the Earth-Moon point for a rendezvous instead of a lunar parking 
orbit. Although this rendezvous scheme has been examined by a number of people since 1960, the first pub- 
lished work was presented by Raithel (Reference 12) in 1966. Because the libration-point rendezvous 
method has several attractive advantages over the lunar orbit rendezvous method, a few comments are in 
order. 


In Figure 8.7, two transfer trajectories between an Earth parking orbit and the L 1 point are 
illustrated. A direct transfer (trajectory #1) uses two impulses, while an indirect transfer (trajectory #2) 
employs a third impulse at the closest approach to the Moon. In both cases, the velocity increment (AF) at 
the Earth parking orbit is approximately equal to the AF that is required for escape velocity. A digital com- 
puter simulation is needed to obtain the remaining AF’s because a patched-conic approximation is not very 
accurate for these calculations. Some computer results of Nicholson (Reference 88) are quoted here. For 
the direct transfer, AF L1 = 2350 fps for a 4-day transit. For the indirect transfer, SAF = 1900 fps (the AF 
near the Moon is about 1000 fps, and AF L1 = 900 fps), but the trip time is 6 days. The same AF's are 
required for the reverse trajectories. 

A transfer between the L 1 point and the Moon is also shown in Figure 8.7 (trajectory #3). In 
this transfer, AF Ll = 1700 fps, and the impact velocity at the Moon is about 7790 fps for a 24-hour transit. 

For landing sites near the lunar equator or for surface stay times of a few days, the total AF 
requirements for lunar orbit rendezvous are slightly lower than the AF requirements for libration-point ren- 
dezvous. However, for landings at higher lunar latitudes or longer stay times, rather large plane changes 
are often necessary, and the AF cost for these plane changes is significantly higher if lunar orbit rendez- 
vous is employed (Reference 89). In some instances, 
this increased AF penalty is severe, and large por- 
tions of the lunar surface are not accessible with the 
present Apollo spacecraft. (An even greater restric- 
tion on lunar surface accessibility is imposed by 
requiring a continuous abort capability for the LLV.) 
On the other hand, the AF’s for the libration-point 
rendezvous technique are virtually identical for any 



Figure 8.7-Transfer paths to and from the L 1 point 
in the Moon's orbital plane. Not to scale. 
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lunar landing site or stay time because the necessary plane changes can be accomplished quite cheaply at 
the Lj point. This low-AF plane-change capability might very well favor a rendezvous at Lj for lunar sur- 

face missions of long duration or high latitude. 

Because the L 1 point is stationary with respect to the lunar surface, the libration-point ren- 
dezvous technique also has an important operational advantage, namely, an infinite launch window for the 
LLV to and from the lunar surface. This timing advantage, along with the low- AT plane-change capability, 
makes the L, point an ideal location for a lunar logistics staging depot (or a lunar-surface rescue facility). 
Supplies could be transported to and from the L 1 point with reusable shuttle vehicles. A low-thrust vehicle 
could be used between the Earth parking orbit and the L 1 point, but a high-thrust vehicle would be needed 
for the transfer between the L 1 point and the lunar surface. 

From the preceding discussion, it is clear that a comprehensive investigation of the libration- 
point rendezvous concept is warranted. This investigation should consider several different transfer modes 
as well as various staging possibilities. The usefulness of the L g point for rendezvous should also be 

examined. 

2. An Interplanetary Transportation System With Terminals at the Sun-Planet Collinear Points 

A modified version of the libration-point rendezvous technique described above may be useful for 
planetary missions. Although a large variety of mission profiles with libration-point rendezvous could be 
contrived, only one possibility is presented here. 

Consider a reusable shuttle vehicle that operates between the L ^ point of the Sun-Earth system 
and the L g point of a Sun-planet system. The outbound transfer path for the shuttle vehicle is shown in 
Figure 8.8; the inbound transfer path is just the mirror image of the outbound transfer path. Either a high- 
or low-thrust rocket could be used to perform the transfer, but it is quite possible that a hybrid system (high- 
and low-thrust) would be more efficient. The transfer of the interplanetary shuttle vehicle (ISV) is initiated 
by applying a small impulse at the Sun-Earth L x point and then starting the low-thrust engine. As the ISV 
passes close to the Earth (point A in Figure 8.8), a much larger impulse is administered. During the helio- 
centric portion of the transfer, only low thrust is employed. At the closest approach to the planet (point B 
in Figure 8.8), another large impulse is applied, and a planetary landing vehicle (PLV) is separated from 
the ISV and lands on the planet.* The ISV then proceeds to the Sun-planet L g point where capture is 
effected by another small impulse. A reverse procedure is used for the inbound transfer. 

Additional shuttle vehicles can be used for the transfer of crew members, fuel, and other supplies 
between the Sun-Earth h x point and the Earth. This phase of the mission is similar to the elliptical orbit 
pickup procedure discussed in Reference 90. 

From the AF standpoint, the libration-point rendezvous technique probably does not have any sig- 
nificant advantage over the elliptical orbit rendezvous method. However, the increased flexibility in the 


*It may be more desirable to separate the PLV earlier or later, depending on the descent mode (atmospheric or 
propulsive braking), guidance requirements, staging ratios, and safety or other factors. 
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Figure 8.8— Transfer path between the Sun-Earth L ^ 
point and a Sun- planet L 2 point. Not to scale. 


timing of various operations* (e.g., rendezvous, 
abort, and landing) may justify a quantitative inves- 
tigation of the libration-point rendezvous concept. 

3. Deep Space Communications Using a Relay 
Satellite at an Earth-Moon Equilateral- 
Triangle Point 

Future deep space missions will require 
much higher rates of data transmission than have 
been attained to date. Laser communications sys- 
tems show considerable promise for furnishing these 
high data rates (Reference 92). However, to avoid 
occupations caused by cloud cover and Earth rota- 
tion, a relay satellite between the deep space vehi- 
cle (DSV) and the Earth-based station will usually 
be required.! A microwave link is used between the 
Earth and the relay satellite, while an optical link 
(laser beam) is used between the relay satellite and 

the DSV. 


The selection of the orbit for the relay sat- 
ellite is very important, and the tradeoffs involved 

in this choice are discussed in Reference 93. Although a preference for a synchronous (24-hour) orbit was 
„„ in Reference 93 the possibility of stationing the relay satellite at an Earth-Moon equilateral- 
triangle libration point was not considered. It is hoped that the equilateral-triangle points will be consid- 
ered in future tradeoff studies because these locations seem to possess certain advantages over a sy 
nous orbit. Some of these advantages are 

(1) Less Earth occultation (almost none) of the optical communications link between the relay 

satellite a^ ^ ^ ugually be excluded from the field of view of the communications receiver on 

the DSV. (The Earth is a strong noise source.) d 

(S) The maximum relative velocity normal to the line of sight between the relay satellite and 

nsv is smaller Therefore, the pointing requirements for the laser beam are less stringent. 

DSV — The cost ' or init . al placement of the relay satellite is reduced by about 0.5 km/sec 

(Reference 15). 

Of course the equilateral-triangle points also have some disadvantages when compared with a synchronous 
orbit (e.g.’, the path loss for the microwave link is greater), and a complete systems evaluation is nee 

before a final choice can be made. 

.The timing problems associated with orbiting or landing missions to the satellites or Jupiter are particularly 
involved (Reference 91). niimhpr of Varth stations, but frequent switchovers would be 

necesZ Ul Btrus:re"iti V o:of a IZuser beam is rather difficult, these switchovers would probably result 
in the loss of some data. 
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The libration-point communications system could also be used for another purpose. With two relay 
satellites, one at L ^ and the other at it would be possible to take advantage of the large separation 
distance (about 6.66 x 10 5 km) to establish an interferometric tracking and navigation system with a long 
baseline (Reference 94). However, because the baseline is rotating with respect to the DSV, it may be ad- 
visable to generate two more baselines by placing a third relay satellite at the Earth-Moon Lg point. (The 
three baselines would form a triangle.) 

B. Scientific Usefulness 

1. A Multiple-Satellite System for Monitoring Solar-Induced 

Phenomena in and Beyond the Earth's Magnetosphere 

Since 1957, scientific satellites have collected a large amount of data concerning particles and 
fields in the Earth's magnetosphere and in the interplanetary medium surrounding it (References 95 to 97). 
These measurements have led to many impressive discoveries of previously unknown phenomena, but future 
progress will be largely dependent on the ability to perform simultaneous observations at different locations. 
Only in this way will it be possible to separate temporal from spatial variations. 

Good coverage of cislunar space has already been achieved by using several satellites in highly 
eccentric Earth orbits. However, beyond the Moon's orbit, the data are usually discrete. To improve this 
situation, Robinson (Reference 14) has suggested that a satellite should be placed in the Earth s orbit about 
6 x 10 6 mi in front of the Earth (this corresponds to a Sun-Earth isosceles-triangle point at about 0=4°) 
where it would continuously monitor the interplanetary medium that will be occupied by the Earth 4 days 
later. Another possibility has been advanced by Meissinger and Greenstadt (Reference 17), who have exam- 
ined the feasibility of obtaining long-term measurements in the Earth's geomagnetic tail by stationing a 
satellite in the vicinity of the Sun-Earth L g point. These stationary satellites would provide some of the 
data that are needed to discern the time variation of solar-induced phenomena in the Earth's neighborhood, 
but more spatial coverage would be desirable. 

Symmetrical spatial coverage could be obtained with a network of four satellites, as shown in Fig- 
ure 8.9. These satellites are stationed at the isosceles-triangle points / 1 and / g (6 = ±2°)* and the collinear 
points L x and L g . t When data collected at these points are compared with measurements taken in cislunar 
space, some interesting correlations may be apparent. It may also be possible to find some correlation with 
Earth weather data. 

Notice that the satellite at the L 1 point could also function as an early warning station for solar 
flares. Although the warning time would only be about 8 to 25 minutes (average velocities of solar flares 
are between 1000 and 3000 km/sec), this may be sufficient for some purposes (e.g., a supersonic transport 
flying at a high altitude could descend to a safer altitude on short notice). 


*The value of 9 that is given here was chosen somewhat arbitrarily. The optimum choice can be determined only 
after considering various tradeoffs (station-keeping cost, communications requirements, scientific value, etc.). 

f The placement of a second satellite at L 2 distance, but along the nominal direction of the geomagnetic tail, 
might also be profitable. By analyzing radio signals between this satellite and the satellite at L 2 , it may be possible 
to monitor the electron de nsit y in a cross section of the geomagnetic tail. The station-keeping cost for the displaced 
satellite would be about |,F cy | = (S L2 “ 1)<V = 1-25 x 10 — & g, with S y = 1.05 x 10 5 km. 
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Figure 8. 9-Multiple- satellite network. Satellites 
anchored at L g , I v and / g . 


The station-keeping costs for the satellites 
at the isosceles triangle points can be obtained from 
Figure 2.8. For 6 = ±2°, the cost is |F C ] = 1.50 x 
10 _6 g. At the collinear points, the costs are 
negligibly small if DSIF tracking is utilized (see 
Chapter V). 

2. Low-Frequency Radio Astronomy from the 
Earth-Moon L g Point 

The low-frequency cutoff for Earth-based 
radio telescopes is about 10 MHz. This limit could 
be extended with an Earth-orbiting radio telescope, 
but it would be extremely difficult to observe fre- 
quencies below 1 MHz because the ionosphere of the 
Earth is a source of low-frequency radio noise. How- 
ever, a radio telescope located at the Earth-Moon L g 
point would always be shielded from the Earth's 
ionosphere by the Moon.* Therefore, the Earth-Moon 


L 0 point may be an ideal location for low-frequency radio astronomy experiments. Of course, possible inter- 
ference from the geomagnetic tail would not be completely eliminated, but the Earth-Moon L g point would 
remain outside of the tail region for periods of almost 3 weeks. 


*The use of the Earth-Moon L z point as a site for radio astronomy was originally suggested by Hornby and Allen 
(Reference 15). 
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CONCLUSIONS AND RECOMMENDATIONS 


The translation control of a satellite in the vicinity of a libration point can be accomplished in a rather 
simple manner. For instar^e, at an unstable collinear point, a single-axis control using only range and 
range-rate measurements is sufficient. Moreover, the station-keeping costs are very low (comparable to 
costs for synchronous satellite station keeping), and are well within the capability of present satellite 
translation control systems. 

The control system can be mechanized in many ways, and some possibilities were presented in this 
study. (The author realizes that the engineering difficulties associated with the cable control technique are 
rather formidable, but the method is intriguing.) The author urges that specific designs of control systems 
for libration-point satellites be the focal point of future work; these studies should include detailed investi- 
gations of propulsion systems and measurement techniques. 

The potential usefulness of a libration-point satellite is apparent (see Chapter VIII). Therefore, the 
early construction and flight testing of a libration-point satellite is strongly recommended. A costly devel- 
opment program may not be necessary because a satisfactory test vehicle could probably be obtained by 
modifying an Applications Technology Satellite. The Earth-Moon L 1 point might be a convenient location 
for station-keeping tests, but a lunar far-side communication capability could also be achieved if the tests 
were conducted in a quasi-periodic orbit around the Earth-Moon L 2 point. 
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APPENDIX A 


DERIVATION OF EQUATION 5.18 


Given the Gaussian density function 
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In G 2 (t), let v = - y 2 /2, then 
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The time average of Equation A7 is 
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where / 0 (u) is a modified Bessel function and u = Ky \op. The time average of G^(t) can be written 
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Integrating by parts, Equation A9 becomes 
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The first portion of Equation A10 vanishes. Making use of the relation 
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the remaining portion of Equation A10 is converted to 
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Therefore, the average control acceleration is given by 
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APPENDIX B 


ECCENTRICITY CORRECTION FOR THE Z-AXIS OSCILLATION 
AT THE EARTH-MOON L 2 POINT 


Following the procedures used in Chapter III, it is readily deduced that the differential equation for the 
eccentricity correction of the z-axis oscillation at a collinear libration point is (see Equation 2.43c) 


where 


Z n2 + B L z n2 ~ 3p B L z n V 
p = - e cos t 


(Bl) 

(B2) 


and 


z ni ~ A zi cos (B3) 

The solution of Equation Bl can be written in the form 

Z n2 ~ A z2 C0S ^z + + A 'z2 C0S ^z ~ (B4) 

At the Earth-Moon Lg point, B L2 = 3.19042, e = 0.05490, and <o z = 1.78618. With these constants, it is 
found that 


A z g - 0.057460 A zl 

and (B5) 

A' zZ = - 0.10416 A zl . 
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